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The  future  development  of  key  areas  in  the  field  of  matrix  structural 
analysis  is  considered.  The  relationship  of  these  areas  to  the  overall 
design  problem  is  established,  and  sources  of  past  difficulties  within  the 
areas  are  identified.  Avenues  of  approach  to  the  difficulties  are  suggested 
and  illustrated  by  reference  to  recently  conducted  research.  In  the  area 
of  modelling,  the  triangular  membrane  element  and  the  bar-panel  ideal¬ 
ization  are  evaluated.  A  method  of  improving  results  obtained  from  the 
analysis  of  models  composed  of  triangular  membrane  elements  is  pre¬ 
sented.  The  loss  of  numerical  accuracy  in  linear  structural  analysis  is 
considered  and  methods  of  improving  accuracy  are  described.  Problems 
in  the  field  of  nonlinear  structural  analysis  are  reviewed.  The  possiblity 
of  an  interactive  computational  approach  to  the  synthesis  of  comples 
structures  is  suggested,  and  an  interactive  graphics  approach  to  the 
synthesis  of  a  structure  subjected  to  steady -state  vibration  is  demonstrated. 
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SECTION  I 
INTRODUCTION 


In  searching  for  a  suitable  theme  for  this  opening  paper,  the  writer  was  impressed  with 
two  features  of  the  conference  agenda  that  seem  especially  important.  The  first  feature  is 
that  the  field  of  interest  is  broad.  We  are  not  only  interested  in  the  ordinary  behavior  of 
structures  subjected  to  loads  well  within  the  elastic  range,  but  we  are  also  concerned  with  the 
exceptional  behavior  resulting  from  the  application  of  loads  outside  this  range.  We  are  interested 
in  the  application  of  matrix  algebra  to  the  analysis  of  structures,  as  well  as  the  process  of 
structural  idealization  that  must  precede  analysis.  We  are  interested  in  structural  synthesis, 
which  is  the  process  of  designing  a  structure  to  meet  a  given  set  of  criteria.  We  are  interested 
in  dynamics  as  well  as  statics,  and  in  the  application  to  our  problems  of  computing  devices, 
such  as  graphics  consoles.  Finally,  no  restrictions  are  imposed  on  the  intended  purposes  of 
the  structures  we  consider.  Many  of  our  methods  can  apply  equally  well  to  aircraft,  missiles, 
space  vehicles,  ships,  submarines,  buildings,  and  bridges. 

A  second  important  feature  of  the  agenda  is  that  the  need  for  progress  in  the  field  is 
great.  Designers  must  cope  with  ever  more  stringent  requirements  resulting  from  projects 
of  an  advanced  nature,  such  as  aircraft  designed  for  vertical  take-off,  supersonic  and  hyper¬ 
sonic  flight,  and  vehicles  intended  for  lunar  missions.  In  striving  for  structural  integrity  and 
efficiency  the  designer  is  not  aided  by  spectacular  breakthroughs,  as  in  other  fields.  Whereas, 
the  speed  of  flight  has  increased  from  about  300  miles  per  hour  to  orbital  velocity  in  30  years, 
a  ratio  of  roughly  fifty  to  one,  the  structural  designer  considers  himself  fortunate  if  the 
strength- to- weight  ratios  of  his  materials  have  increased  by  50  percent  in  the  same  period  of 
time. 

Because  the  field  is  so  broad,  only  certain  key  areas — where  the  need  for  accelerated 
progress  seems  most  urgent — have  been  selected  for  discussion.  The  objective  of  this  paper 
is  to  consider  the  future  development  of  these  key  areas  in  order  to: 

•  establish  the  significance  of  each  area  as  related  to  the  overall  problem  of  structural 
design 

•  identify  clearly  the  sources  of  past  difficulties 

•  suggest  some  approaches  to  these  difficulties 

•  illustrate  certain  approaches  by  reference  to  research  recently  conducted  by  the  writer. 
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Engineers  and  scientists  working  in  a  highly  technical  field,  such  as  structural  analysis, 
often  concentrate  their  efforts  in  specialized  areas.  It  is  hoped  that  this  paper  will  serve  to 
provide  a  useful  overview  of  a  significant  fraction  of  the  entire  field,  and  that  it  will  point  the 
way  to  avenues  of  research  that  will  prove  useful  in  the  future.  These  areas  include  discrete 
element  idealization,  linear  structural  analysis,  nonlinear  analysis,  synthesis  and  the  inter¬ 
active  application  of  remote  graphics  consoles  to  problems  of  steady- state  vibration. 


SECTION  n 

DISCRETE  ELEMENT  IDEALIZATIONS 


This  section  is  limited  to  a  consideration  of  discrete  idealizations,  or  models,  composed 
of  finite  numbers  of  elements.  These  elements  are  connected  at  joints  or  nodes  which  also 
serve  as  application  points  for  idealized  discrete  equivalents  of  actual  external  loads. 

Two  principal  classes  of  models  are  defined.  A  displacement  model  is  composed  of  ele¬ 
ments  having  force-deformation  characteristics  such  that  continuity  of  displacements  exists 
across  element  boundaries,  while  the  characteristics  of  equilibrium  models  are  such  that 
continuity  of  stresses  is  maintained  (Reference  1). 

The  present  status  ofmodellingcanbe  summarized  roughly  as  follows:  Framed  structures 
can  be  modelled  without  difficulty.  Structural  elements  that  represent  plates,  shells  and  solid 
bodies  have  been  devised.  These  elements  provide  reliable  results  only  if  properly  employed 
by  experienced  personnel.  Better  plate  and  shell  elements,  together  with  better  methods  of 
interpreting  results  derived  from  these  elements,  are  needed.  Better  methods  of  evaluating 
the  usefulness  of  elements  are  also  desired. 

METHODS  OF  ANALYSIS 

A  clear  distinction  should  be  made  between  the  structural  model  and  the  method  whereby 
the  model  is  analyzed.  Such  a  distinction  is  necessary  to  provide  a  proper  basis  for  discussing 
the  difficulties  involved  in  the  two  areas. 
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The  static  behavior  of  a  linear  discrete  structure  is  governed  by  equations  of  equilibrium 
and  compatibility.  Three  methods  of  solving  these  equations  are  defined.  In  the  Force  Method, 
displacements  are  eliminated  from  the  governing  equations.  The  resulting  equations  are  solved 
for  element  forces,  which  serve  as  the  basis  for  calculating  displacements.  In  the  Displacement 
Method,  element  forces  are  eliminated  from  the  governing  equations.  The  resulting  equations 
are  solved  for  displacements.  Elements  forces  are  calculated  from  these  displacements. 
In  the  Unified  Method,  the  equations  are  solved  in  no  predetermined  order. 

The  Force  and  Displacement  Methods  are  subsets  of  the  Unified  Method  as  demonstrated 
by  Kosko  (Reference  2).  The  definitions  of  these  methods  do  not  depend  on  the  structural 
model.  Any  of  the  three  methods,  therefore,  can  be  employed  to  analyze  either  displacement 
or  equilibrium  models.  The  results  obtained  through  the  use  of  any  of  the  three  methods  are 
the  same,  except  for  loss  of  accuracy  during  computation,  when  applied  to  any  particular  model 
regardless  of  type,  since  the  governing  equations  being  solved  in  each  case  are  the  same. 

The  degree  of  accuracy  inherent  in  a  given  model  should  not  be  attributed  to  a  particular 
analysis  method.  For  example,  the  statement  that  the  Force  Method  gives  upper  bounds  for 
deflection  influence  coefficients  is,  generally,  not  true, 

MODELLING  CRITERIA 

An  attempt  to  assess  the  usefulness  of  structural  models  must  be  based  on  criteria  of 
some  kind.  Following  are  three  criteria  that  seem  to  be  generally  accepted: 

1.  Convergence  -  As  the  sizes  of  elements  in  a  model  approach  zero,  the  errors  in 
computed  stresses  and  deflections  should  approach  zero.  The  convergence  should 
be  monotonic. 

2.  Rate  of  Convergence  -  The  usefulness  of  an  element  is  proportional  to  its  rate  of 
convergence. 

3.  Computational  Efficiency  -  The  usefulness  of  an  element  is  proportional  to  the  ratio 
of  accuracy  to  computational  effort. 

INTERPRETATION 

The  stresses  and  deflections  mentioned  in  the  criteria  are  nominal  values  derived  from 
the  model  according  to  some  rule  of  interpretation.  For  example,  nodal  deflections  are 
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normally  Interpreted  to  be  structural  deflections.  Various  options  are  open  for  stress  inter¬ 
pretation,  Thus,  the  stresses  at  the  center  of  a  membrane  element  can  be  interpreted  as  the 
stresses  at  the  correspondingpoint  of  the  structure.  Some  interpretation  rule  must  be  establish¬ 
ed  before  the  modelling  criteria  can  be  applied. 

MODEL  EVALUATION 

One  method  for  evaluating  discrete  membrane  elements  is  to  consider  a  plane  stress 
problem  having  a  known  solution.  For  example,  aplate  of  a  given  shape  is  subjected  to  boundary 
stresses.  A  closed  form  solution  for  the  resulting  internal  stress  distribution  is  known.  The 
plate  is  then  idealized  with,  say,  triangular  elements  of  a  certain  size.  Stresses  and  deflections 
are  numerically  computed  and  compared  with  the  known  solution.  The  process  is  repeated  for 
elements  of  different  sizes.  Convergence  is  studied. 

A  METHOD  OF  DISCRETE  ELEMENT  EVALUATION 

The  following  paragraphs  describe  a  different  approach  to  element  evaluation.  It  is 
similar  to  the  preceding  approach,  except  that  after  the  plate  is  modelled  with  discrete  ele¬ 
ments  it  is  analyzed  mathematically  instead  of  numerically.  Difference  equations  are  written 
for  stresses  and  deflections.  These  equations  are  solved  and  the  results  are  compared  with 
known  solutions  to  provide  expressions  for  stress  and  deflection  errors  as  functions  of  panel 
size.  This  method  has  several  advantages.  First,  a  precise  expression  for  rate  of  convergence 
is  obtained.  Conclusions  about  convergence  are  not  obscured  by  errors  that  occur  when  stresses 
and  deflections  are  numerically  computed.  Second,  results  are  obtained  for  a  continuous  range 
of  element  sizes.  Third,  a  study  ofthe  error  expressions  can  lead  to  a  better  understanding  of 
element  behavior.  Asa  consequence,  methods  of  correcting  the  results  or  of  improving  the 
elements  themselves  can  be  devised.  Fourth,  the  best  rules  for  interpreting  results  can  be 
established. 
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SECTION  m 

THE  EQUILATERAL  TRIANGULAR  CONSTANT  STRESS  MEMBRANE 
ELEMENT,  POISSON’S  RATIO  =  1/3 


Figure  1  shows  a  triangular  element  applicable  to  the  Matrix  Force  method.  The  forces 
F1,  F2>  and  Fg  act  parallel  to  the  element  edges.  The  corresponding  deformations  are  edge 
elongations.  Stresses  in  the  element  are  assumed  to  be  uniform.  Appendix  1  describes  the 
element  in  more  detail. 


Figure  1.  Triangular  Membrane  Element 


LINEAR  STRESS  DISTRIBUTION,  TRIANGULAR  ELEMENTS 


Appendix  1  presents  a  study  of  the  behavior  of  the  triangular  element  in  a  general  linear 
stress  field.  The  element  is  considered  to  be  part  of  a  grid  of  the  kind  shown  in  Figure  2, 
although  the  grid  shape  may  vary.  The  grid  is  a  model  of  an  isotropic  plate  of  uniform  thick¬ 
ness.  Loads  act  upon  the  boundaries  distributed  in  such  a  manner  that  internal  stresses  are 
linear  functions  of  the  coordinates. 
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With  these  assumptions  the  analysis  of  the  model  is  a  simple  matter.  The  true  displace¬ 
ments  of  the  plate  are  first  calculated.  These  displacements  are  then  applied  to  joints  of  the 
model  and  the  resulting  external  loads  are  determined.  Element  forces  and  element  stresses 
are  computed.  Appendix  1  contains  the  details  of  the  analysis.  The  results  follow.  These  results 
also  apply  to  the  uniformly  stressed  triangular  element  employed  in  the  Matrix  Displacement 
method. 

Edge  Loads 

The  analysis  shows  that  normal  and  tangential  loads  applied  to  a  joint  on  the  edge  of  the 
idealized  plate  should  be  equal  to  the  boundary  normal  and  shear  stresses,  respectively, 
acting  at  the  joint  multiplied  by  the  length  of  the  edge  of  the  element  and  the  thinkness  of  the 
plate.  No  loads  are  applied  to  internal  j oints.  The  idealized  plate  is  in  equilibrium  under  these 
loads  when  the  stress  distribution  is  linear. 

Deflections 

The  computed  deflections  are  exact  when  the  loads  are  applied  as  described  in  the  pre¬ 
ceding  paragraph.  Deflections  of  joints  of  the  idealized  structure  should  be  interpreted  as 
deflections  at  corresponding  points  of  the  plate. 

Stresses 

An  attempt  was  made  to  find  a  point  in  each  element  of  the  model  at  which  the  computed 
uniform  stresses  in  the  element  are  equal  to  the  true  stresses  at  the  corresponding  point  in 
the  plate.  No  such  point  exists.  Next,  an  attempt  was  made  to  find  a  point  in  each  element  at 
which  the  sum  of  the  squares  of  the  errors  in  the  stresses  is  minimized.  Discovery  of  such 
a  point  would  provide  a  rule  for  stress  interpretation.  No  single  point  was  found.  Instead  the 
location  of  the  optimum  point  was  found  to  depend  on  the  type  of  linear  stress  field  acting  on 
the  plate. 

Since  no  unique  optimum  point  was  found,  a  study  was  made  of  errors  at  the  centroid  of 
the  triangle.  The  results  are  given  by  Equation  (33)  and  are  summarized  as  follows:  First, 
the  stress  errors  at  the  centroid  of  the  element  converge  to  zero  as  the  length  of  the  side  of 
the  element  approaches  zero.  Second,  the  errors  are  proportional  to  the  length  of  the  side  of 
the  element  and  the  stress  distribution  parameters  b^,  b b^,  and  b^.  Third,  the  errors  are 
independent  of  element  location. 
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These  results  indicate  that  calculations  based  on  the  triangular  element  generally  con¬ 
verge,  as  long  as  the  stress  field  contains  no  singularities,  because  the  field  in  the  region  of 
an  element  will  tend  to  approach  linearity  as  the  element  becomes  smaller. 

Percentage  Error 


The  error  in  a  quantity  measured  or  calculated  by  some  approximate  procedure  can  be 
expressed  as  a  percentage  of  the  known  correct  value.  This  procedure  is  not  applicable  for 
present  purposes  because  small  errors  in  stresses  that  have  very  small  or  zero  exact  values 
produce  large  percentage  errors  which  are  misleading.  The  following  procedure  is  adopted: 
Errors  in  tensile  stress  components  are  divided  by  an  equivalent  tensile  stress,  ~cr  ,  defined 
as  the  uniform  tensile  stress  that  would  produce  the  same  amount  of  strain  energy  in  the  plate 
as  the  actual  stress  distribution.  Errors  in  shear  stress  components  are  likewise  divided  by 
an  equivalent  uniform  shear  stress,  r  .defined  in  a  similar  manner.  These  reference  stresses, 
which  are  constant  for  any  particular  loading  condition,  are  derived  as  follows:  The  strain 
energy  stored  in  a  plate  subjected  to  plane  stress  is: 

W:T  [°x  +  °y  ~ZvaxCTy  +2<l  +  *')TxZy  ]  dA  (I) 


where  t  =  plate  thickness,  E  =  Young’s  modulus,  v  =  Poisson’s  ratio  and  cr  ,  cr  ,  and  r 

_  x*  y’  xy 

are  components  of  the  stress  tensor.  If  cr  -  cr  (constant)  and  cr  =  r  =0,  then 

x  y  xy 


C T 


If  r  =  r  (constant)  and  cr  =  cr  =0,  then 
xy  x  y 


r  = 


EW 


A(  1  +v  1  J 2  ( I  +  v  ) 

Stress  Correction  Procedure  for  Constant  Stress  Triangular  Elements 


(2) 


The  notion  of  the  stress  field  in  the  region  of  an  element  becoming  approximately  linear 
as  the  element  decreases  in  size  suggests  a  method  of  improving  stresses  computed  from 
triangular  elements.  The  method  is  based  on  the  assumption  that  a  linear  stress  field  that 
approximates  the  computed  stress  distribution  can  be  found  in  the  region  of  an  element. 
Appendix  2  contains  the  derivation  of  the  method. 


When  the  correction  procedure  is  applied  to  the  approximate  stresses  derived  from  the 
equilateral  triangular  grid  subjected  to  a  linear  stress  field,  exact  results  are  obtained.  The 
method  also  produces  excellent  results  in  the  following  highly  nonlinear  case. 

A  concentrated  load  acts  upon  one  vertex  of  the  triangular  plate  shown  in  Figure  3. 
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Figure  3.  Loaded  Triangular  Plate 

Equilibrium  is  produced  by  stresses  acting  along  the  lower  edge.  The  following  stresses  and 
displacements  of  the  plate  are  derived  from  Reference  4. 
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The  plate  was  replaced  by  structural  models  composed  of  equilateral  triangular  membrane 
elements,  and  stresses  computed  from  the  models  were  compared  with  exact  results  and  with 
stresses  computed  according  to  the  correction  procedure  of  Appendix  2.  Three  grid  sizes  were 
studied:  fine,  medium,  and  coarse.  Figure  4  shows  the  medium  fineness  model.  The  fine 

grid  had  twice  as  many  elements  per  edge  and  the  coarse  grid  had  half  as  many.  The  triangular 


Figure  4.  Medium  Grid  Model 


plate  was  truncated  as  shown  because  the  concentrated  load  applied  at  the  vertex  produces 
infinite  stresses  and  deflections,  at  that  point,  which  are  not  realistic. 

Each  model  was  analyzed  in  the  following  manner.  Nodal  displacements  were  accurately 
calculated  from  Equation  4.  The  matrix  of  displacements  was  multiplied  by  the  assembled 
stiffness  matrix  to  produce  external  loads.  Element  deformations  were  calculated  from  dis¬ 
placements,  and  element  stress  components  (  orx  oy  rXy)  were  computed  from  element 
deformations.  The  model  was  thus  analyzed  as  a  kinematically  determinate  structure.  This 
procedure  has  three  advantages  compared  to  idealizing  edge  loads  and  computing  deflections 
and  stresses  from  these  loads.  First,  the  loss  of  accuracy  involved  in  solving  simultaneous 
equations  is  avoided.  Second,  the  uncertainty  resulting  from  idealizing  edge  loads  according 
to  some  arbitrary  rule  is  avoided.  Third,  a  study  of  the  computed  external  loads  may  suggest 
rules  for  load  idealization.  The  kinematically  determinate  structure,  although  a  somewhat 
unusual  subject  for  study,  should  yield  reliable  information  on  convergence. 
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Figure  5  shows  the  convergence  of  external  load  at  the  node  indicated  by  a  circle  in 
Figure  4.  The  error  in  the  load  acting  at  the  node,  expressed  as  a  percentage  of  the  load  acting 


p 

% 


0  12  3  4 

Length  of  Element  Edge 
(  I nch e  s  ) 

Figure  5.  Convergence  of  Interior  Joint  Loads 

at  the  vertex  of  the  triangular  plate,  is  plotted  as  a  function  of  length  of  element  edge.  The 
load  should  be  zero  since  the  node  is  on  the  interior  of  the  plate.  Figure  5  indicates  a 
satisfactory  convergence. 

Figure  6  shows  the  convergence  of  shear  stress  at  the  point  indicated  by  a  plus  mark  {+)  in 
Figure  4.  The  shear  stress  error,  as  a  percentage  of  the  mean  shear  stress  in  the  plate,  is 
plotted  as  a  function  of  length  of  element  edge.  The  corrected  shear  stress  shows  a  very 
satisfactory  convergence. 

The  shear  and  normal  loads  acting  at  nodes  on  the  upper  and  lower  edges  of  the  plate  were 
found  to  approach  the  shear  and  normal  stresses  acting  at  these  nodes  multiplied  by  the  length 
of  the  element  edge  and  the  thickness  of  the  plate. 

The  results  shown  in  Figures  5  and  6  are  typical  of  the  entire  plate.  The  errors  in  the 
medium-grid -corrected  stresses  are  on  the  order  of  three  percent  compared  to  about  20  per¬ 
cent  for  uncorrected  stresses.  The  corresponding  values  for  the  fine  grid  are  about  one  percent 
and  nine  percent. 
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Figure  6.  Convergence  of  Shear  Stress 

The  equilateral  triangular  element,  on  the  basis  of  convergence  of  external  loads  and 
corrected  stresses,  can  be  considered  to  meet  the  modelling  criteria  in  a  satisfactory  manner. 
A  good  deal  of  confidence  can  be  placed  in  the  convergence  of  deflections  when  loads  are  speci¬ 
fied,  as  a  result  of  the  satisfactory  convergence  of  loads  when  deflections  are  specified.  This 
conclusion  holds,  provided  that  excessive  accuracy  is  not  lost  in  the  solution  of  simultaneous 
equations. 
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SECTION  IV 

THE  RECTANGULAR  BAR-PANEL  IDEALIZATION 

Figure  7  shows  a  bar-panel  idealization  of  a  plate.  The  model  consists  of  axially  loaded 
bars,  and  plate  elements  or  panels  that  carry  shear  only.  Panels  are  connected  to  bars  at 


Figure  7.  Bar- Panel  Model 


midpoints  of  panel  edges.  The  point  of  connection  between  a  panel  and  a  bar  is  referred  to  in 
the  following  as  a  shear  node.  The  point  of  intersection  between  two  bars  is  called  a  joint. 
Normal  edge  loads  act  on  joints,  and  tangential  edge  loads  act  on  shear  nodes.  Bars  need  not 
be  equally  spaced.  Bars  and  panels  are  constant  stress  elements. 


LINEAR  STRESS  DISTRIBUTION,  BAR-PANEL  MODEL 

The  behavior  of  the  bar-panel  idealization  in  a  linear  stress  field  is  analyzed  in  Appendix 
3,  accounting  for  Poisson’s  ratio  and  considering  unequally  spaced  bars,  with  the  following 
results: 


Edge  Loads 

As  in  the  case  of  triangular  elements,  the  normal  load  applied  to  a  joint  should  be  equal 
to  the  boundary  normal  stress  acting  at  the  joint,  multiplied  by  the  distance  between  midpoints 
of  adjacent  panels  and  the  plate  thickness.  The  tangential  load  applied  to  the  shear  node  should 
be  equal  to  the  boundary  tangential  stress  at  the  node,  multiplied  by  the  length  of  the  panel 
edge  and  the  plate  thickness. 
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Deflections 

Deflections  at  joints  are  exact  for  the  linear  stress  field  when  loads  are  applied  as  in  the 
preceding  paragraph.  Deflections  of  shear  nodes  should  be  disregarded. 

Stresses 

Normal  stresses  in  bars  should  be  interpreted  as  normal  stresses  at  points  corresponding 
to  joints  in  the  plate.  Shear  stresses  in  panels  should  be  interpreted  as  shear  stresses  at 
points  corresponding  to  panel  centers  in  the  plate.  Interpreted  in  this  manner,  the  stresses 
produced  by  the  bar-panel  idealization  are  exact  for  the  linear  stress  field. 

QUADRATIC  STRESS  DISTRIBUTION,  BAR-PANEL  MODEL 

The  study  of  the  behavior  of  the  bar-panel  model  in  the  linear  stress  field  throws  no  light 
on  rate  of  convergence  because  stress  and  deflection  errors  are  zero.  The  analysis  of  model 
behavior  in  a  general  quadratic  stress  field  is  therefore  described  in  Appendix  3.  Because  of 
the  complexity  of  this  analysis,  two  simplifying  assumptions  are  made.  Bar  elements  are 
assumed  to  be  equally  spaced,  and  Poisson’s  ratio  is  assumed  to  be  zero. 

The  stress  function  for  the  quadratic  distribution  involves  four  independent  parameters 
d^,  dg,  d^,  and  dg.  The  effects  of  the  distributions  resulting  from  parameters  d^  and  d4  on 
the  bar-panel  model  are  essentially  the  same  as  the  effects  of  the  distributions  resulting  from 
d^  and  d^,  because  of  isotropy  of  the  model.  Only  the  distributions  d^  and  d^  need  be  considered. 
The  results  are  as  follows: 

Edge  Loads 

For  thed1  distribution,  edge  loads  should  be  idealized  in  the  same  manner  as  for  the  linear 
stress  field.  For  the  d^  distribution  this  rule  produces  edge  loads  that  do  not  satisfy  equilibrium 
of  moments.  In  Appendix  3  this  moment  imbalance  is  eliminated  by  a  balancing  couple  uniformly 
distributed  along  the  edges  of  the  plate  parallel  to  the  Y-axis. 

Deflections 

The  errors  in  joint  deflections  are  proportional  to  the  squares  of  lengths  of  panel  edges. 
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Stresses 


Errors  in  stresses  for  the  distribution  are  zero.  For  the  d2  distribution,  errors  in 
tensile  stresses  are  zero.  The  error  in  the  shear  stress  is  proportional  to  the  square  of  the 
length  of  the  panel  edge  parallel  to  the  x-axis. 

In  order  to  derive  quantitative  results  for  the  rate  of  convergence,  consider  the  square 
plate  shown  in  Figure  8  subjected  to  the  quadratic  stress  field  defined  by  putting  d0  =  1,  d1  - 
d4  =  d&  =  0.  From  Equation  8 5  the  stresses  are  0-^  =  0,  cry=  6xy,  =  -3x  .  From  Equations 
2  and  85  the  reference  shear  stress  is 


From  Equations  6  and  102  the  absolute  value  of  the  percentage  error  is 


1  r  m 

where  m  is  the  number  of  shear  elements  per  side  of  the  plate.  Figure  9  shows  the  convergence 
as  a  function  of  m. 


y 


Figure  8.  Square  Plate 
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Number  of  elements  per  side 

Figure  9.  Convergence  of  Shear  Stress  in  Bar-Panel  Model 

The  bar-penel  idealization  can  be  considered  to  meet  the  modelling  criteria  in  a  satis¬ 
factory  manner  not  only  on  the  basis  of  the  preceding  analysis  but  also  on  the  basis  of  many 
past  comparisons  with  test  results  (Reference  11). 
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SECTION  V 

LINEAR  STRUCTURES 


The  element  force  and  displacement  responses  of  structures  in  this  category  are  linear 
functions  of  loads,  support  displacements,  and  unassembled  element  deformations. 

BASIC  APPROACHES  TO  STATIC  ANALYSIS 

The  Force,  Displacement,  and  Combined  Methods  differ  in  the  order  in  which  unknowns 
are  computed,  but  all  three  methods  provide  solutions  to  the  same  governing  equations.  All 
three  methods  are  applicable  to  any  linear  discrete  model.  The  essential  difference  is  in  order 
of  solution.  This  difference  at  first  may  not  appear  to  be  significant,  but  disputes  concerning 
the  comparative  suitability  of  the  three  methods  in  various  applications  have  persisted.  A  brief 
review  of  the  pros  and  cons  of  this  matter  is  offered. 

Required  Unknowns 

If  only  deflections  are  required,  the  Displacement  Method  has  an  advantage  because  the 
effort  needed  to  compute  forces  can  be  avoided.  Similarly  if  only  forces  are  required,  the 
Force  Method  has  an  advantage. 

Structural  Model 

The  Displacement  Method  has  an  advantage  when  the  model  has  a  relatively  low  ratio  of 
number  of  nodal  displacements  to  number  of  redundants.  The  Force  Method  has  an  advantage 
when  the  ratio  is  high.  The  first  situation  tends  to  occur  when  the  structure  has  many  members 
in  parallel,  the  second  when  the  structure  has  many  members  in  series.  Particular  types  of 
elements  seem  better  suited  to  one  method  than  the  other.  For  example,  the  rectangular  con¬ 
stant  stress  membrane  element  seems  better  suited  to  the  Displacement  Method  than  the  bar- 
panel  idealization  because  the  latter  requires  additional  nodes  at  points  where  panels  are  joined 
to  bars.  Lumping  (joining  elements  that  serve  parallel  functions)  is  a  useful  procedure  in  the 
Force  method  for  reducing  the  number  of  redundants. 

Conditioning  of  the  Governing  Equations 

Large  complex  structural  models  have  a  tendency  to  yield  governing  equations  that  are 
difficult  to  solve  without  significant  loss  of  accuracy.  The  order  of  solution  can  be  decisive 
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in  successfully  solving  such  equations.  Forsyth  has  discussed  this  matter  in  some  detail 
(Reference  5).  As  an  example,  consider  the  following  matrix  equation: 


A  solution  is  sought  on  the  basis  of  a  Jordanian  elimination  process.  The  most  straight¬ 
forward  approach  is  to  attempt  to  eliminate  the  unknowns  in  the  order  x^,  Xg,  x^,  x^  to  yield 

a  solution  for  x_,  using  diagonal  elements  as  pivots.  The  following  equation  results  after 
0 

elimination  of  x^,  Xg  and  x^: 


"-I00000I  1000000 

1000000  -100000 1 


(9) 


This  equation  is  poorly  conditioned.  If  a  solution  is  attempted  with  a  computer  that  carries 
less  than  seven  significant  decimal  digits,  the  method  fails. 


If  the  unknowns  are  eliminated  in  the  order  x^,  x^,  x^,  x^,  Xg,  with  elements  (3,5),  (1,4), 
(4,1),  (5,3)  and  (2,2)  of  the  coefficient  matrix  as  pivots,  then  the  following  results  can  be 
obtained  with  the  aid  of  a  two  digit  machine. 

x  x  x  it  1  :  (  0.50  0.50  -0.50  0.50  0.50  )  (10) 

l  I  2  3  4  5  J  l  > 

This  result  is  a  very  good  approximation  to  the  exact  solution  which  is 

{*,  ««  s  «,  *,}  =  nn? '■  -|H0"6'  '■  l+,°6 1  (l" 

The  scalar  equations  contained  in  Equation  8  are  actually  well-conditioned  in  the  sense  that 
an  accurate  solution  can  be  obtained  as  a  result  of  carrying  a  small  number  of  significant 
digits,  provided  the  order  of  solution  is  optimized.  It  seems  worthwhile  to  distinguish  between 
such  equations  and  equations  that  are  poorly  conditioned.  The  term  “order  sensitive”  is  pro¬ 
posed  to  denote  simultaneous  equations  having  the  property  that  solution  accuracy  is  sensitive 
to  order  of  elimination. 
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The  example  shows  that  some  sets  of  equations  are  order  sensitive.  An  approximate 
optimum  order  of  solution  can  be  discovered  through  the  use  of  the  procedure  known  as  pivot 
selection  (References  5  and  6),  A  suitable  pivot  selection  strategy  incorporated  into  equation 
solving  routines  for  matrix  structural  analysis  can  significantly  improve  the  accuracy  of 
computed  results. 


The  influence  of  order  of  elimination  on  solution  accuracy  has  a  bearing  on  the  question 
of  Force  versus  Displacement  and  Unified  Methods.  In  the  Unified  Method  no  predeter¬ 
mined  restrictions  are  placed  on  order  of  solution.  The  Unified  Method  in  conjunction  with 
pivot  selection  therefore  provides  in  every  case,  a  higher  probability  of  success  than  either 
of  the  other  two  methods.  Indications  are  that  the  governing  equations  are  usually  well 
conditioned,  although  order  sensitive,  if  the  structure  is  stable.  The  Unified  Method  with 
pivot  selection  should  be  a  very  reliable  procedure.  The  principal  disadvantage  of  the  method 
is  that  all  of  the  governing  equations  must  be  processed  in  one  batch. 

The  other  two  methods  involve  prescribed  orders  of  solution.  No  one  has  been  able  to 
show  which  order  has  the  greatest  potential  for  accuracy.  However,  it  is  known  that  a  structural 
model  with  one  or  more  very  stiff  members  produces  poorly  conditioned  Displacement  Method 
equations  unless  certain  nodal  displacements  are  deleted,  whereas,  the  Force  Method  equations 
are  well-conditioned  in  such  cases.  The  detrimental  effect  of  stiff  members  plus  the  lack  of 
suitable  pivot  selectionproceduresmaybethecauseof  difficulties  sometimes  experienced  with 
Displacement  Method  solutions. 

The  Force  Method  involves  the  difficulty  of  redundant  selection.  One  approach  to  this 
problem  is  to  base  the  selection  on  experience  and  intuition.  This  approach  is  virtually  useless 
for  most  large  structures.  A  second  approach  is  to  select  redundants  by  comparison  with  a 
standarized  model,  or  a  model  for  which  a  mathematical  analysis  is  available.  This  approach 
should  be  reliable  when  a  suitable  comparison  model  can  be  found.  Many  structures  are  too 
complex  to  be  dealt  with  in  this  manner.  A  third  approach  to  redundant  selection  is  known  as 
the  Structure  Cutter.  This  approach  is  based  upon  the  application  of  pivot  selection  procedures 
to  solve  the  equilibrium  equations  in  terms  of  the  redundants  plus  the  application  of  column 
weighting  factors  proportional  to  member  stiffnesses  (Reference  6).  This  procedure  eliminates 
the  redundant  selection  problem  by  completely  automating  the  selection  process.  The  Structure 
Cutter  has  the  additional  advantage  of  optimizing  the  order  of  solution  of  the  equilibrium  and 
compatibility  equations  involved  in  the  Force  Method. 
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Experience  with  both  the  Force  and  Displacement  Methods  has  led  the  writer  to  believe 
that  the  Force  Method  with  optimized  redundant  selection  tends  to  provide  a  greater  degree  of 
accuracy  than  the  Displacement  Method. 

Computational  Efficiency  and  Convenience 

The  influences  upon  computational  efficiency  of  the  configuration  of  the  structural  model 
and  the  particular  unknowns  to  be  calculated  have  been  discussed.  Computational  efficiency 
is  also  dependent  upon  the  design  of  the  computer  program  and  the  skill  of  the  programmer. 
Extremely  efficient  programs  applicable  to  special  classes  of  statically  indeterminate 
structures  have  been  devised.  User  convenience  is  a  function  of  the  degree  of  automation  built 
into  the  program.  Maximizing  the  degree  of  automation  saves  valuable  engineering  time  and 
increases  reliability.  These  comments  apply  equally  to  the  Force  and  Displacement  Methods. 

PROGRAMMING  SYSTEMS 

A  valid  approach  to  developing  an  automated  system  of  structural  analysis  is  to  write 
programs  applicable  to  particular  types  of  structures  such  as  beams,  plates,  and  shells  of 
revolution.  A  second  approach  is  to  develop  general  purpose  programs.  This  approach  has 
proven  feasible  and  worthwhile  because  many  components  of  typical  aircraft,  missile,  and 
space  vehicle  structures  do  not  fall  into  recognizable  categories.  The  following  paragraphs 
describe  a  general  purpose  system  recently  developed  at  the  Douglas  Aircraft  Company  under 
sponsorship  of  the  Air  Force  Flight  Dynamics  Laboratory, 

FORMAT  System 

The  FORTRAN  Matrix  Abstraction  Technique  (FORMAT)  is  a  general  purpose  program¬ 
ming  system  for  matrix  structural  analysis.  Detailed  documentation  of  FORMAT  is  contained 
in  References  7,  8,  and  9.  A  description  of  some  applications  of  the  system  is  given  in 
Reference  10.  The  purpose  here  is  to  describe  briefly  the  project  objectives,  the  programming 
concepts,  and  the  degree  to  which  the  objectives  have  been  met. 

The  major  criteria  originally  identified  for  FORMAT  were  generality,  communicability, 
flexibility,  universality,  and  efficiency  (Reference  6).  Generality  meant  that  the  program  should 
be  applicable  to  a  broad  range  of  structures.  Communicability  meant  that  the  mathematical 
procedures  should  be  formulated  in  such  a  manner  as  to  be  intelligible  to  engineers  working 
in  related  fields.  Flexibility  meant  that  users  who  were  expert  in  matrix  structural  analysis 
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would  be  able  to  alter  procedures  to  meet  changing  requirements  with  no  additional  programming 
effort.  Universality  meant  machine  independence.  Efficiency  meant  maximum  output  at  minimum 
engineering  and  computing  cost. 

Some  of  these  criteria  tend  to  be  contradictory  and  compromises  are  necessary.  An 
effort  was  made  to  optimize  the  overall  effectiveness  of  the  program. 

The  FORMAT  system  is  programmed  in  three  phases.  Phase  I  generates  matrices  from 
raw  data.  The  output  of  Phase  I  includes  matrices  required  for  the  following  functions:  analysis 
of  structures  subjected  to  thermomechanical  loading,  analysis  of  symmetric  structures, 
element  modification,  substructure  joining,  and  analysis  of  elastic  stability  and  undamped 
vibration.  The  phase  also  generates  joint  coordinates  of  some  commonly  used  structural 
shapes,  and  maintains  and  updates  files  of  data. 

Phase  II  performs  matrix  operations  according  to  abstraction  instructions  written  by 
the  user.  A  typical  instruction  is  FXO  =  FX.  MULT.  XO,  which  means  that  matrix  FXO  is 
equal  to  FX  postmultiplied  by  XO.  The  user  thus  has  freedom  in  determining  the  sequence  of 
operations.  The  program  processes  abstraction  instructions,  prints  diagnostics,  allocates 
matrices  to  available  tapes  in  an  optimized  manner,  and  executes  instructions.  Matrices  up 
to  order  2000  can  be  processed. 

Phase  III  prints  Force  and  Displacement  Method  reports  and  provides  a  graphical  output 
capability. 

The  original  criteria  for  the  system  have  essentially  been  met.  The  Displacement  and 
Force  method  generators  provide  a  selection  of  displacement  and  equilibrium  elements  suf¬ 
ficient  to  represent  a  broad  range  of  engineering  structures.  The  matrix  notation  and  extensive 
user  and  programming  documentation  provide  communicability.  The  matrix  abstraction 
language  provides  flexibility.  The  machine  independence  provided  by  the  FORTRAN  coding  is 
such  that  the  program  has  been  successfully  run  on  a  number  of  machine  configurations, 
including  the  IBM  7094,  IBM  7044-7094  direct  couple,  UNIVAC  1108,  GE  635,  and  IBM  360/65. 
Some  degree  of  computing  efficiency  was  sacrificed  to  meet  the  other  criteria,  but  the  savings 
in  programming  effort  and  maintenance  resulting  from  generality,  flexibility,  and  relative 
machine  independence  are  considered  to  more  than  offset  possible  increases  in  computing  cost. 
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SECTION  VI 

NONLINEAR  STRUCTURES 


The  governing  equations  for  structures  in  this  category  are  nonlinear  in  the  unknown  dis¬ 
placements  and  element  forces.  The  nonlinearities  result  from  any  form  of  behavior  that  does 
not  accord  with  the  assumptions  of  small  deflections  and  linear  elasticity,  hi  the  presence  of 
creep,  responses  to  static  loads  are  time-dependent.  The  search  for  solutions  to  nonlinear 
structural  problems  involves  certain  difficulties. 

BASIC  LAWS 

The  first  difficulty  is  that  the  basic  laws  of  material  behavior  are  not  completely  known. 
The  strain  components  of  an  incremental  volume  of  material  are  functions  of  stress  components, 
past  histories  of  these  components,  and  histories  of  other  environmental  factors  such  as 
temperature.  The  laws  of  macroscopic  material  behavior  are  probably  different  for  each  alloy. 

REPETITIVE  SOLUTIONS 

The  quadratic,  cubic,  and  biquadratic  equations  in  a  single  variable  are  the  only  nonlinear 
equations  having  closed  form  solutions.  All  other  nonlinear  equations,  including  the  governing 
simultaneous  equations  for  nonlinear  structures,  must  be  solved  by  repetitive  processes. 
These  processes  are  usually  iterative  or  stepwise  applications  of  linear  approximations. 

In  the  iterative  approach,  a  solution  is  sought  for  a  particular  set  of  values  of  the  independent 
variables.  Various  structural  parameters,  such  as  element  stiffnesses,  are  functions  of  depen¬ 
dent  variables,  such  as  stress,  which  are  unknown.  Values  of  the  parameters  are  assumed. 
Approximate  dependent  variables  are  calculated  and  revised  parameters  are  determined  on  the 
basis  of  the  new  independent  variables.  The  process  is  repeated  until  values  of  the  unknowns 
are  obtained  that  satisfy  all  of  the  governing  equations  to  the  required  degree  of  accuracy. 

In  the  stepwise  approach  the  response  is  assumed  to  be  linear  over  small  increments  of 
the  independent  variable,  which  can  be  load  or  time.  The  structural  parameters  existing  at  the 
beginning  of  the  interval  are  assumed  to  apply  to  the  entire  interval.  Various  refinements 
of  the  technique  exist.  The  iterative  and  stepwise  approaches  can  be  combined. 
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In  both  the  iterative  and  stepwise  approaches  convergence  is  important.  The  results 
must  converge  to  fixed  values  as  the  number  of  iterations  increases,  or  as  the  size  of  steps  in 
the  stepwise  method  decreases.  Convergence  can  be  deceptive.  The  results  can  appear  to 
converge  without  being  even  approximately  correct.  Theprocesscan  be  very  slowly  convergent 
or  divergent.  The  results  of  the  stepwise  process  can  have  properties  which  are  characteristics 
of  the  method  of  solution  rather  than  of  the  structure. 

MULTIPLE  ROOTS 

The  possession  of  multiple  roots  is  a  property  of  nonlinear  equations.  For  a  given  set 
of  loads  several  different  sets  of  structural  responses  can  satisfy  the  governing  equations. 
Only  one  response  is  appropriate  to  the  prior  history  of  loading.  The  iterative  approach  con¬ 
verges  on  roots  that  depend  upon  the  choice  of  the  initial  approximation  and  the  method  of 
solution.  The  stepwise  approach  has  the  advantage  of  reproducing  the  actual  loading  history 
and  consequently  has  the  potential  capability  of  always  converging  on  the  correct  set  of  roots. 

NONVALIDITY  OF  SUPERPOSITION 


A  third  difficulty  is  that  the  responses  of  a  nonlinear  structure  to  two  separate  loading 
conditions  cannot  be  added  to  produce  the  response  to  the  sum  of  the  two  conditions.  The 
solution  for  each  condition  must  be  separately  computed.  The  unit  load  approach  as  a  labor 
saving  device  for  developing  responses  to  actual  loading  conditions  is  not  applicable. 

COMPUTING  EFFORT 

A  fourth  difficulty  is  that  the  magnitude  of  the  computing  task  is  very  significantly  in¬ 
creased  for  nonlinear  structures  as  a  result  of  the  required  repetitive  approach  and  the 
inapplicability  of  superposition. 

APPROACHES  TO  THE  DIFFICULTIES 

One  possible  means  of  remedying  the  lack  of  information  on  material  behavior  is  to  con¬ 
tinue  study  of  the  subject  on  the  molecular  and  crystalline  levels.  The  forces  that  govern  the 
behavior  of  these  bodies  may  be  simpler  to  classify  and  describe  than  the  behavior  of  macro¬ 
scopic  bodies  such  as  tensile  coupons. 
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A  second  approach  is  to  continue  deliberate  efforts  to  establish  categories  of  idealized 
material  behavior  that  are  useful  in  structural  analysis.  Categories  already  established  in¬ 
clude  nonlinear  elasticity,  perfect  platicity  and  various  time -dependent  stress-strain  laws. 
The  criterion  here  is  not  how  closely  natural  behavior  is  approximated,  but  the  usefulness  of 
the  approximation. 

The  need  for  repetitive  solutions,  the  existence  of  multiple  roots  and  the  inapplicability 
of  superposition  are  fundamental  and  unavoidable.  The  impact  of  these  difficulties  can  be 
diminished  by  searching  for  optimum  methods  of  solution  that  minimize  computing  effort  and 
maximize  the  probability  of  converging  on  the  correct  set  of  roots. 

The  irreducible  minimum  amount  of  calculation  unfortunately  can  be  expected  to  be  at 
least  an  order  of  magnitude  greater  than  the  effort  required  for  the  linear  case,  which  already 
has  a  reputation  as  a  formidable  computing  task.  A  requirement  therefore  exists  for  vastly 
improved  computing  capability.  Perhaps  this  improvement  can  be  realized  by  developing 
specialized  computers  for  the  solution  of  nonlinear  simultaneous  equations.  The  requirements 
in  the  field  of  structural  analysis  may  be  sufficient  to  justify  such  development. 

Finally,  a  need  exists  to  define  carefully  the  objectives  of  nonlinear  structural  analysis 
in  view  of  the  difficulties  involved.  The  basic  objective  for  airframe  and  space  vehicle  design 
is  to  find  methods  that  are  useful  for  design  purposes,  not  to  model  nature  in  great  detail. 
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SECTION  vn 

STRUCTURAL  SYNTHESIS 


Structural  synthesis  is  the  procedure  of  designing  a  structure  to  meet  a  criterion  of 
merit  within  constraints.  For  example,  a  structure  can  be  designed  to  have  minimum  weight 
(the  criterion  of  merit),  and  to  have  the  capability  of  supporting  multiple  loading  conditions 
in  such  a  manner  that  stresses  and  deflections  do  not  exceed  allowable  values  (the  constraints). 

Minimum  weight  is  often  employed  as  a  criterion  of  merit.  Maximum  cost  effectiveness  is 
usually  a  more  realistic  criterion  although  more  difficult  to  apply. 

Some  examples  of  constraints  are  as  follows: 

•  No  structural  element  shall  fail  by  fracture  or  buckling  at  ultimate  load 

•  No  yielding  shall  occur  at  limit  load 

•  No  general  instability  shall  occur 

•  Deflections  must  not  exceed  allowable  values 

•  Aeroelastic  effects  (e.g.,  loss  of  aileron  effectiveness)  shall  not  exceed  allowables 

•  No  flutter  shall  occur 

•  No  fatigue  failure  shall  occur  during  the  design  life. 

A  list  two  or  three  times  as  long  can  easily  be  prepared.  These  constraints  are  important 
design  requirements  that  must  be  satisfied. 

The  independent  variables  of  the  problem  include  loads  and  thermal  environment.  The 
unknown  dependent  variables  include  structural  concept,  structural  geometry,  member 
sizes,  materials,  means  of  attachment,  and  thermal  protection. 

NONLINEARITY  OF  THE  SYNTHESIS  PROBLEM 

Figure  10  shows  a  simple  pin-jointed  statically  indeterminate  truss.  The  problem  is  to 
determine  the  member  cross-sectional  areas  A^,  A2»andAg  in  such  a  manner  as  to  minimize 
weight,  subject  to  the  constraints  that  no  member  shall  have  a  negative  area.  The  loading 
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conditions  are  as  follows:  in  condition  1,  (p ^  =  l«<jb 2  =  ™  conditio11  2,  =  0,  <£>2  =  1.  All 

members  have  the  same  allowable  tensile  stress  cr  The  stress-strain  relationships  are 

a 

linear  to  failure.  The  material  density  is  unity. 


Analysis  of  the  problem  yields  the  following  relationships, 
for  Condition  1: 

' 

<x0  (20.48  A2  A3  4  10  A3  A,  h  10  A,  A2)>  (30  A2  4  30  A3)  ( 

aQ  (20.48  Az  A3  +10  A3  A,  +  10  A,  Ag  )>  38.4  A3<£, 

crQ  (2C  48  Ag  A3  +  10  A3  A,  4  10  A,  A2)  >  38.4  A2  c 

for  condition  2: 


(12) 


<r  (20.48  Az  A  4  10  A3  A,  +  10  A(  A2)>(40A3  -  40  A,,  )</>,,  > 

cr  (20.48  A2  A3  h  10  A3  A,  +  10  A,  A2)>  (51.2  A3+  50  A,  )<^)2 
Q  \  I  o 

<ra  (20.48  A2  A3  +  10  A3  A  (  +-  10  A,  Aa»>(5i.2  .Ag  4-  50  A,  )  <p2 
Also  A(>0  A2>  O. 

The  mathematical  formulation  is  in  terms  of  inequalities  rather  than  equations.  These 
inequalities  are  quadratic  in  the  unknowns  A^,  A2  and  A^,  even  though  the  structure  is 
linear.  As  more  structural  elements  are  added,  the  formulation  of  the  synthesis  problem 
becomes  more  highly  nonlinear  (not  merely  quadratic)  and  extremely  complex. 
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NATURE  OF  SYNTHESIS 


The  general  synthesis  problem,  which  is  the  problem  designers  must  solve,  is  formidable 
when  stated  in  mathematical  terms.  The  following  complexities  are  involved:  First,  realistic 
criteria  of  merit,  such  as  maximum  cost  effectiveness,  are  not  easy  to  define  with  precision. 
Second,  the  list  of  constraints  is  very  long.  Third,  the  number  of  variables  to  be  considered  is 
large.  Fourth,  the  problem  is  nonlinear.  Fifth,  the  mathematical  formulation  tends  to  produce 
inequalities  rather  than  equations.  Finally,  the  synthesis  problem  includes  all  analytic 
problems,  such  as  linear  stress  and  deflection  analysis,  plasticity,  buckling,  and  vibration. 
The  mathematical  formulation  of  the  general  synthesis  problem  is  many  orders  of  magnitude 
more  complex  than  the  most  difficult  problem  in  the  field  of  structural  analysis.  The  possibility 
of  producing  a  general  computer  program  for  structural  synthesis  seems  extremely  remote. 

APPROACH  TO  SYNTHESIS 


Subsets  of  the  general  problem  have  been  studied.  Working  computer  programs  have  been 
developed  (References  12  and  13)  These  programs  are  useful  tools  and  they  are  significant 
advances  in  the  state  of  the  art.  They  contribute  to  the  solution  of  the  general  problem  but  they 
do  not  solve  it,  yet  this  problem  must  be  and  is  being  solved  daily  to  some  degree  of  approxi¬ 
mation,  Two  questions  suggest  themselves:  First,  howis  the  problem  currently  being  handled? 
Second,  how  can  the  capability  of  the  computer  and  the  methods  of  matrix  analysis  be  more 
efficiently  applied  to  structural  synthesis? 

Obviously  the  answer  to  the  first  question  is  that  the  designer  does  the  best  he  can  with 
the  tools  he  has.  These  tools  include  digital  computers  and  computer  programs,  but  the 
designer’s  most  valuable  aid  is  the  ability  to  think  creatively  about  the  problem.  The  computer 
is  entirely  lacking  in  this  ability  which  the  human  brain  has  in  abundance.  Even  before  the 
invention  of  the  slide  rule  men  were  able  to  do  a  creditable  job  of  structural  design  aided 
solely  by  computations  performed  with  pencil  on  paper.  On  the  other  hand,  the  computer 
is  far  superior  to  man  in  the  area  of  complex  arithmetic  calculation. 

Recognition  of  this  complementary  relationship  between  man  and  machine  led  to  the  con¬ 
cept  of  Computer  Aided  Design  (CAD)  (References  14  and  15).  According  to  this  concept  man 
and  computer  are  placed  in  close  contact,  so  that  interaction  time  is  minimized  and  maximum 
benefits  are  gained  from  the  combination  of  abilities.  The  concept  is  implemented  through  the 
use  of  remote  consoles,  especially  graphic  consoles,  connected  to  large  computers  provided 
with  operating  systems  that  allow  the  user  almost  instant  response. 
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The  idea  of  applying  CADprinciples  to  synthesis  of  complex  structures  is  worth  consider¬ 
ation  in  view  of  the  difficulties  of  developing  a  completely  automated  system.  Unfortunately 
many  computing  tasks  associated  with  matrix  structural  analysis  are  too  large  for  the 
interactive  approach.  On  the  other  hand,  some  tasks  can  be  approached  in  this  manner.  Such 
an  application  is  described  in  the  next  section. 


SECTION  VIII 

AN  INTERACTIVE  APPROACH  TO  STEADY  STATE  VIBRATION 

Voluminous  output  tends  to  be  characteristic  of  matrix  structural  analysis.  Various  means 
such  as  tables  of  maxima  and  minima  are  employed  to  reduce  the  amount  of  data  the  user  must 
examine.  Computer  graphics  provides  an  extremely  effective  means  of  condensing  output.  In  a 
recent  dynamic  response  study  (16)  an  interactive  graphics  console  displayed  a  moving  picture 
of  a  vibrating  structure.  Such  an  application  is  especially  useful  because  the  dynamic  display 
conveys  an  impression  of  mode  shapes  that  cannot  be  easily  gained  from  inspecting  printed 
data. 


The  following  paragraphs  describe  another  application  of  the  graphics  console  to  structural 
vibration. 

In  this  application,  the  steady  state  responses  of  a  damped  structure  to  oscillating  inputs 
were  studied.  As  a  result  the  matrix  formulations  presented  in  Appendix  4  were  developed. 
These  formulations,  which  are  generalizations  of  the  Displacement,  Force,  and  Unified 
Methods  of  static  analysis,  yield  rigorous  solutions  to  the  equations  of  motion.  The  formulations 
are  applicable  to  linear  discrete  structures.  Each  structural  element  has  its  own  elastic  and 
damping  constants.  Masses  are  attached  to  structural  joints.  All  kinds  of  support  and  zero 
support  are  considered.  Forcing  functions  that  are  sinusoidal  functions  of  time  of  a  single 
frequency  but  differing  phases  are  applied.  These  forcing  functions  are  external  loads,  support 
displacements,  and  unassembled  element  deformations.  The  unassembled  deformations  provide 
the  basis  for  rapid  computation  of  structural  changes.  Equations  are  given  for  time-dependent 
force  and  displacement  responses  of  the  structure. 
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In  static  analysis  the  effects  of  changing  a  structural  element  can  be  simulated  by  intro¬ 
ducing  into  the  element  fictitious  unassembled  deformations.  The  magnitudes  of  these  defor¬ 
mations  are  such  that  the  assembled  forces  and  deformations  of  the  modified  element  are  the 
same  as  they  would  be  if  the  element  stiffness  actually  were  changed.  This  idea  has  been 
adapted  to  the  analysis  of  vibrating  structures.  The  resulting  procedure  for  calculating 
modified  response  involves  much  less  effort  than  computing  the  effects  of  a  change  from  basic 
data. 


The  computations  take  place  in  two  steps.  The  first  step  takes  place  in  the  batch  mode. 
In  this  step  the  response  of  the  original  structure  is  calculated.  In  addition  certain  matrices 
are  computed  which  are  required  in  the  second  step  as  a  basis  for  computing  the  effects  of 
structural  changes. 

In  the  second  step  the  user  inputs  through  the  graphics  console  the  identification  of  the 
elastic  element,  mass,  or  damper  he  wishes  to  modify.  He  also  inputs  the  factor  by  which  he 
wishes  to  multiply  the  stiffness,  the  mass,  or  the  damping  constant.  Changes  to  the  deflection 
response  are  calculated  and  displayed.  The  user  can  modify  the  elements  in  any  order  by  any 
amount.  The  same  element  can  be  modified  repeatedly. 

Figure  11  shows  a  structure  consisting  of  a  horizontal  gridwork  of  beams  and  a  number 
of  vertical  axially  loaded  members,  as  displayed  on  the  face  of  the  cathode  ray  tube.  A  mass 
is  attached  to  each  node,  and  an  oscillating  force  of  a  given  frequency  acts  at  the  upper  end  of 
the  center  vertical  bar.  Each  element  has  its  own  damping  constant. 

Figure  12  shows  the  vibrating  structure  at  a  particular  instant  of  time.  A  moving  picture 
also  can  be  displayed.  The  display  can  be  rotated  about  each  of  the  coordinate  axes,  and  slices 
of  the  structure  can  be  displayed  for  improved  visibility. 
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Figure  12.  Response  of  a  Damped  Elastic  Structure 

Figure  13  shows  the  list  of  elements  that  are  available  for  modification.  The  tables  show 
the  original  value  and  the  latest  value  of  the  modified  element  parameters.  When  the  user 
requests  a  change  to  a  structural  element,  1-1/2  minutes  are  required  for  the  machine  to 
compute  and  display  the  revised  deflection  pattern.  The  modification  program  is  thus  a  true 
interactive  application  of  the  graphics  console  to  matrix  structural  analysis. 
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ELEMENT 

TYPE 

SEQ. 

NO. 

ELEM. 

NO. 

ORIGINAL 

VALUE 

SCALING 

FACTOR 

LAST  VALUE 

ELASTIC 

1 

1 

.  17500E  04 

. 10000E  01 

. 17500E  04 

ELASTIC 

2 

2 

.  17500E  04 

. 10000E  01 

. 17500E  04 

ELASTIC 

3 

3 

. 17500E  04 

. 10000E  01 

. 17500E  04 

ELASTIC 

4 

4 

.  17500E  04 

. 10000E  01 

.17500E  04 

ELASTIC 

5 

5 

. 17500E  04 

. 10000E  01 

.  17500E  04 

MASS 

6 

87 

.33333E  00 

. 10000E  01 

. 33333E  00 

MASS 

7 

88 

. 33333E  00 

. 10000E  01 

. 33333E  00 

MASS 

8 

89 

.33333E  00 

. 10000E  01 

. 33333E  00 

MASS 

9 

90 

.33333E  00 

. 10000E  01 

. 33333E  00 

MASS 

10 

91 

. 33333E  00 

. 10000E  01 

. 33333E  00 

DAMPING 

11 

92 

. 17500E  04 

. 10000E  01 

. 61250E  04 

DAMPING 

12 

93 

. 61250E  04 

.10000E  01 

. 61250E  04 

DAMPING 

13 

94 

. 61250E  04 

, 10000E  01 

. 61250E  04 

DAMPING 

14 

95 

. 61250E  04 

. 10000E  01 

. 61250E  04 

DAMPING 

15 

96 

. 61250E  04 

. 10000E  01 

.  61250E  04 

SEQ. 

ORIGINAL 

LAST  VALUE 

SCALING 

NEW  VALUE 

NO. 

VALUE 

FACTOR 

1 

. 17500E  04 

. 17500E  04 

100.0 

Figure  13.  Modified  Element  Parameters 
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SECTION  DC 
SUMMARY 


The  field  of  structural  modeling  can  benefit  from  the  establishment  of  generally  accepted 
criteria  for  discrete  elements,  and  systematic  methods  of  element  evaluation.  The  evaluation 
methods  should  be  as  free  as  possible  from  the  influences  of  extraneous  factors,  such 
as  loss  of  accuracy  in  the  process  of  numerical  evaluation.  An  analytic  evaluation  avoids  this 
difficulty  and  provides  a  better  understanding  of  element  behavior. 

Loss  of  accuracy  in  solving  the  governing  structural  equations  has  always  been  a  problem. 
The  accuracy  obtained  in  solving  the  equations  by  an  elimination  process  can  be  sensitive  to 
the  order  of  elimination  even  when  the  equations  are  well  conditioned.  Accuracy  can  be 
improved  by  exercising  care  in  establishing  the  elimination  order  and  by  employing  a  pivot 
selection  technique  during  elimination. 

The  analysis  of  nonlinear  structures  involves  the  difficulties  of  repetitive  solution, 
multiple  roots,  and  nonvalidity  of  superposition.  These  difficulties  significantly  increase  the 
scope  of  the  computing  task.  Methods  that  minimize  this  task  are  needed.  Development  of 
specialized  computing  devices  for  solving  nonlinear  simultaneous  equations  may  be  justified 
for  this  and  other  applications. 

The  scope  of  the  computingtask  involved  in  structural  synthesis  can  be  very  large  because 
a  single  synthesis  can  require  the  solution  of  many  different  kinds  of  analytic  problems,  such 
as  stress  analysis,  stability,  vibration,  and  fatigue.  The  task  can  be  performed  most  efficiently 
through  a  combination  of  the  creative  powers  of  the  analyst  and  the  arithmetic  and  logical 
capability  of  the  computer.  A  suitable  combination  can  be  effected  by  means  of  interactive 
graphics  terminals  and  other  devices  that  bring  the  user  into  close  contact  with  the  computer. 
Such  an  approach  to  the  synthesis  of  complex  structures  should  be  considered.  An  interactive 
application  of  the  remote  graphics  console  to  damped  steady  state  structural  vibration  was 
demonstrated. 
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APPENDIX  I 

EVALUATION  OF  THE  EQUILATERAL  TRIANGULAR  CONSTANT  STRESS 
ELEMENT,  POISSON’S  RATIO  =  1/3 


Figure  14.  Equilateral  Triangular  Element 


The  flexibility  matrix  for  the  triangular  element  shown  in  Figure  14(a)  is 

4 

4 

yi 

4 

/3  . 

where  E  is  Young’s  modulus  and  t  is  the  plate  thickness. 


Figure  14(b)  shows  a  triangular  frame  composed  of  axially  loaded  bars  of  cross-sectional 
area  -/T at/4,  where  a  is  the  length  of  the  edge  of  the  element. 

The  flexibility  matrix  for  the  frame  is  also  given  by  Equation  14,  consequently  the 
equilateral  triangle  can  be  replaced  by  the  framework  of  axially  loaded  bars.  This  substitution 
simplifies  the  discussion. 
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Difference  Equations  for  a  Model  Composed  of  Equilateral  Triangular  Elements 


x 


Figure  15.  Displaced  Bar 

Figure  15  shows  a  bar  of  the  triangular  frame.  The  symbols  for  horizontal  and  vertical 
displacements  of  the  ends  of  the  bar  are  as  shown.  The  subscript  i  indicates  that  the  bar  is  the 
ith  element  of  the  grid.  The  deformation  of  the  bar  is 

e.  =  (u.  -u  )  cos  ctj  +  (v.  -v  )  sin  a.  (15) 

i  i  i  ai 

where  cr  is  the  angle  between  the  bar  and  the  X  axis. 

The  element  force  in  the  bar  is  given  in  terms  of  the  element  deformation  by 

73 

F i  =  4“  Ete i  U6) 


50 


A  FFDL-TR- 68-150 


Figure  16.  Joint  Equilibrium 


Figure  16  shows  an  interior  joint  of  the  idealized  plate,  with  the  adjoining  triangular 
elements.  The  external  loads  Px  and  are  shown  acting  on  the  joint.  Also  acting  on  the 
joint  are  reactions  to  the  element  forces  shown  acting  on  the  corners  of  the  elements  adjacent 
to  the  joint.  The  joint  is  in  equilibirum. 


Px  =  -  Z  F.  COS  Otj 
i 

py  =  "I  Fj  sin  a- 


(17) 


where  the  summations  include  all  element  forces  acting  on  the  joint.  Similar  equations  apply 
to  joints  on  the  edges  and  corners  of  the  plate. 


Equivalent  Uniform  Stresses 

The  flexibility  matrix  Equation  14  is  based  on  the  assumption  of  uniform  stress  in  the 
panel.  The  equivalent  uniform  stresses  are  given  by 


i 

4/3'  /3  /3 

°x 

ot 

3  3  3 

Fi 

0  Jz  Jz 

Fz 

T 

0  -1  1 

F, 

xy  _ 

3 

(18) 


51 


AFFDL-TR-68-150 


This  transformation  is  applicable  when  the  edge  of  the  element  (or  the  bar  of  the  equivalent 
frame)  upon  which  the  F^  force  acts  is  parallel  to  the  x  axis. 

Linear  Stress  Field 


If  the  plate  is  subjected  to  constant  stress,  a  model  composed  of  constant  stress  elements 
yields  exact  results.  A  stress  distribution  one  degree  more  complicated  is  a  distribution  in 
which  stresses  are  linear  functions  of  coordinates.  Consider 


cx  =  2b3  x  +  6b4  y+  2c 


3  \ 


°y  =  6bi  *  +"  2b2  y  +■  2c( 


r  -  -2b„  x  -  2  b  y  -  c 
xy  2  3  2 


(19) 


These  stresses  are  derived  from  the  stress  function 

<p  =  b,  x3  +  b2x2  y+  b3  xy2  +•  b4  y3  +  c,  xZ  +  cg  xy  +  c3  y2  (20) 

The  distribution  given  by  19  is  the  most  general  linear  distribution  that  satisifies 
equilibrium  and  the  condition  V4<£  =0.  The  components  of  the  displacement  of  a  point  having 
the  coordinates  x  and  y  are  given  by 

Eu  =  -  b(  (3vxZ  3y2  )  -  b2(2j/xy )  +b3(x2  -  2y2  -  i/y2  )  +  b^  (6xy) 

-  c  ,(2i/x  )  -  c2<  I  +  v  )y  +  c3(2x)  I 

\  (21) 

Ev=b  (6xy)  +  b  (y2  -2xZ-i/xZ)”b  (2i/xy)-b  (3j/yZ+3x2)  I 

|  2  O  / 

+  c(  (2y  )  -  c2(l+J/)  x  -c3(2i/y  ) 

Two  orientations  of  the  triangular  element  are  considered  as  shown  in  Figures  17  and  18. 
The  figures  also  show  the  notation  for  the  coordinates  of  the  vertices  of  the  elements.  The 
type  II  orientation  ca  n  be  derived  from  the  Type  I  orientationby  changingthe  sign  of  “a”. 
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Substituting  the  coordinates  of  the  vertices  of  the  triangles  shown  in  these  figures  into 
Equations  21,  and  setting  y  -  1/3  yields  expressions  for  the  displacements  of  the  vertices. 
Substituting  these  displacements  into  Equation  15  gives  the  following  expressions  for  the 
elongations  of  the  sides  of  the  triangle  in  the  Type  I  orientation. 


i 


E  -  b  r(-2  X|  -a  )  ■+  b2  (-- 1-  y,  )  +•  b3  (2x(  +  a  )  +  b4  (6y  ) ,  +c,  (-  y)  +  c3(2) 
E-f  =  t4*!  +  3°H-b2(-5^  +  y,+a) 

(4) -3  W 

E4i  =  "l  (4"l  +  0  1  +  b2  (-^3^-  xl  +4  ».)  +  b3  -«) 

+«, 


(22) 
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Expressions  for  the  Type  II  orientation  can  be  obtained  by  replacing  x1  by  x^,  y  by 
y2  and  a  by  -a. 

External  Loads 


Eliminating  F.  among  Equations  16  and  17  gives 

Px  = - El  I  e  .  cos  a. 

I 

73"  _ 

P„  =  -  —~~z  Et  V  e.  sin  a. 

y  4  4-i  1 


(23) 


These  expressions  can  be  evaluated  for  an  interior  joint  such  as  the  one  shown  in  Figure 
16  by  substituting  the  appropriate  element  deformations  given  by  Equations  22  or  similar 
expressions  applicable  to  the  Type  II  orientation. 


Px  =0  | 

Py  =  0  ) 


Figure  19.  Edge  Joint  Equilibrium 
For  an  edge  joint  of  the  kind  shown  in  Figure  19  the  results  are 

Px  =  ot  (2b2xQ  -F2b3  y0  +  c2  ) 

► 

Py  =  «»t(-6b|  x0-Zb2y0-2C|  ) 

where  the  coordinates  of  the  joint  are  x  ,  y  . 

0  0 


(25) 
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It  can  be  seen  from  Equations  25  that  P^  and  Py  are  equal  to  the  shear  and  normal  stresses 
(with  appropriate  signs)  at  the  joint  multiplied  by  the  length  of  the  element  edge  and  the  plate 
thickness. 

This  method  of  establishing  idealized  normal  loads  on  edge  joints  appears  to  neglect 
couples  resulting  from  the  fact  that  normal  edge  stresses  are  linearly  distributed  rather  than 
constant.  Study  shows  that  these  couples  are  exactly  balanced  by  couples  acting  on  other  edges 
that  are  similarly  neglected. 

Stresses 


Substituting  Expressions  22  for  e^  e2>  and  eg  into  Equation  16  yields  expressions  for  the 
element  forces  F^,  Fg,  and  Fg.  Substituting  these  expressions  into  Equation  18  gives  the 
following  expressions  for  the  nominal  uniform  element  stresses: 

°x  =  b2  (“^  °)  +  b3  (2xi  4  b4  (6*i  )  +  c3  (  2  ) 

^  =  b,  (6x,  +-  3a)  +  b2(2y,  f/3c)  4  c,  (2)  ,  {26) 

fxy  'b,  (-T“)  +  V-2V“ltb!(-!yl  - 


Subtracting  the  nominal  stresses  from  the  true  stresses  given  by  Equations  19  and  substituting 
s  -  (x  -  x^)/a  and  t  =  (y  -  y^) /a  leads  to  the  following  matrix  equation  for  stress  errors: 


Aax 

= 

2  b  3 

1  6h4 

1 

s  ” 

4 

r  .  -l 

3  b2  b  3 

Aay 

6b  [  2b 

1  |  2 
t 

t 

-  3b  (  ~-fz  b2 

i - 

> 

X 

1 _ 

CM 

-Q 

1  J 

*r  b<  +  b2  +vf-  b3 

(27) 


where  Act  =  cr  -  cr  and  Acr  =  cr  -  <r  .  Let 

xxx  y  y  y 


A  cr  = 

i 

X 

b 

;  a  = 

_2b3 

j  6b<T 

;  s  = 

's " 

;  b  = 

■/z  Y  .  ~i 

3  b2  b  3 

b* 

<3 

6bl 

|  2b2 

t 

-3b,  -Jz  b2 

At 

xy 

~”2b  2 

!"2b3 

■Jz  73 

IT  b.  +  b2  4  g"  b3 

(28) 
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4‘Ao'  =  AS  +  B  (29) 

Note  that  s  and  tare  nondimens  ional  coordinates  relative  to  the  element  vertex  whose 
coordinates  are  x^,  y^.  To  find  a  point  in  the  element  where  the  nominal  element  stresses  are 
equal  to  the  true  plate  stresses,  set  Act  =  0. 

•'•AS  =  -  B  (30) 

Unfortunately  this  matrix  equation  represents  three  scalar  equations  in  only  two  unknowns. 

2  2  2 

Therefore  no  such  point  can  be  found.  A  point  such  that  Act  +  Act  +  At  is  minimized 

x  y  xy 

can  be  found  in  the  following  manner: 

AT  AS  =  -  AT  B  (31) 

T  - 1  T 

S  -  —  (  A  A  )  A  B  (32) 

Equation  32  was  solved  for  four  sets  of  values  of  the  stress  parameters.  The  results  are 
summarized  in  Table  1. 


TABLE  1 


Case 

Stress  Parameters 

Coordinates  of 
Optimum  Point 

Stress  Component 
Errors 

■9 

b2 

b3 

b4 

s 

t 

Acrx 

EM 

Arxy 

a 

■ 

0 

0 

0 

-J— - 

t 

0 

0 

mm 

b 

1 

0 

0 

c 

1 

2 

2 

3 

0 

■ 

c 

0 

1 

0 

1 

2 

S 

4 

0 

0 

H 

d 

0 

0 

0 

1 

s 

0 

0 

0 

Stress  Errors  at  the  Panel  Centroid 


The  nondimens  ional  coordinates  of  the  panel  centroid  are  s  =  1/2,  t  =  -fW/d>.  Therefore 
from  Equation  29, 


i 

> 

o 

~  ^  b4 

Act 

=  a 

273“ 

D  x 

y0 

3  2 

1 - 

c> 

1 

7E  b  +  7l  b 

2  i  6  3 

where  AcrXo,  Acr^  and  ^Txy0  are  stress  errors  at  the  element  centroid.  These  errors  are 
proportional  to  the  length  of  side  of  the  element  and  independent  of  the  location  of  the  element. 
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APPENDIX  II 

LINEAR  STRESS  CORRECTION  PROCEDURE  FOR  CONSTANT 
STRESS  MEMBRANE  ELEMENTS 


Figure  20.  Plate  Element 

Consider  a  plate  subjected  to  plane  stress.  Replace  the  plate  by  a  model  consisting  of 
constant  stress  elements.  Figure  20  shows  one  such  element,  not  necessarily  triangular. 
Element  forces  are  related  to  element  deformations  by 

F  =  D"1  •  {34) 

where  0  is  the  flexibility  matrix  for  the  element. 


Coordinate  Axes  are  passed  through  point  0  as  shown.  The  components  of  the  stress  tensor, 

<t  „,  cr  and  ~  (constant  in  the  element),  are  obtained  by  a  transformation  of  the  element 
71  y  xy 

forces. 

<r  =  TF  (35) 


where 


and  T  is  a  transformation  matrix. 


(36) 


57 


AFFDL-TR- 68-150 


Assume  that  the  idealized  plate  has  been  analyzed,  so  that  cr  is  known.  The  nominal 

stresses  cr  ,  cr  and  "F  are  approximations  to  the  components  of  the  true  stress  tensor 
^  y  xy 

at  0.  Improved  approximations  are  desired.  To  obtain  such  an  improvement: 

1.  Assume  that  the  stress  distributions  in  the  plate  in  the  region  of  the  element  is 
approximately  linear,  so  that  the  components  of  the  stress  tensor  are  given  by  Equations  19. 
The  corresponding  displacements  in  the  region  of  the  element  are  given  by  Equations  21. 

2.  Calculate  the  deformations  of  the  element  resulting  from  the  displacements  of  Equations 
21  in  the  form 

e  =  eB  B  +  »c  C  (37) 


where 


B 


i 

b  2 
b. 


C  = 


(38) 


and  e  u  and  e  are  coefficient  matrices  whose  elements  are  functions  of  element  geometry. 

-D  U 


3.  Calculate  the  element  forces  from  Equation  34  and  the  nominal  stresses  from  Equation  35. 


=  °B  B''~trC  C 

where 

O’  =  TD  a 

B  B 

and 

Fc  =  TD"  *C 

Let  erQ  be  a  matrix  of  the  true  stresses  at  point  0,  i.e.. 


(39) 

(40) 


(41) 


(42) 
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From  Equation  19, 


(43) 


If  the  stresses  are  constant  (  B  =  0  ),  the  nominal  stresses  in  the  element  given  by  the 
matrix  <r  must  be  equal  to  the  true  stresses  given  by  the  matrix  <tq  since  the  element  is 
a  constant  stress  element.  Therefore,  from  Equation  39 

(44) 

(45) 


and 


°-0 


°C  C 


°B  B 


Equation  45  provides  a  means  of  correcting  the  nominal  stresses  obtained  from  the  use  of 
a  constant  stress  element  since  it  gives  the  true  stresses  at  the  panel  center  in  terms  of  the 
nominal  stresses.  The  equation  is  exact  for  a  linear  stress  distribution.  For  a  nonlinear 
distribution  the  equation  will  provide  an  approximate  means  of  correcting  the  nominal  stresses. 
In  the  case  of  a  nonlinear  distribution  the  matrix  B  can  be  considered  to  be  a  matrix  of 
coefficients  of  a  linear  stress  field  which  approximates  the  nonlinear  field  in  the  region  of 
the  element.  The  matrix  (from  Equation  40)  is  a  known  function  of  element  geometry  and 

stiffness. 


The  matrix  B  can  be  established  by  determining  the  coefficients  b.  so  that  the  stress 
distribution  of  Equations  19  approximately  reproduces  the  stresses  which  have  been  computed 
in  adjacent  elements  at  points  1,  2  ...  n  (Figure  20). 


The  true  stresses  atpoint  i  can  be  equated  to  the  linearly  distributed  stresses  of  Equations  19. 


0  ”i  0  [Ex.  i  6  y . 

i  ,  .  ,  . 

b  | 

~2c3" 

°Vi 

6x,  [  2y.  10  |  0 

1  1  f 

b2 

+ 

2c, 

rxy; 

1  1  1 

0  |  -2x.  l-2y.  |  0 

b3 

~c2 

L  J 

_  1  1  1 

b4 

(46) 


£ T.  :  OK.  B  +  ff 

i  t  o 


1*1,2...  n 


(47) 
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and  a  is  a  scalar  introduced  to  nondimens ionalize  the  coordinates.  The  n  matrix  equations 
given  by  Equations  47  can  be  written 


a  KB  =  cr'  -  cr0'  (51) 

Equation  51  represents  3n  scalar  equations  in  the  unknowns  b^,  bg  and  b^.  The  value 
of  n  should  be  at  least  equal  to  two,  so  that  Equation  51  involves  more  scalar  equations  than 
unknowns.  A  linear  stress  distribution  which  provides  a  least  squares  fit  to  the  true  stresses 
at  points  0,  1,  2,  ...  n  can  be  obtained  as  follows: 

aKTKB  =  KT  (a‘  -  <r'0  )  (52) 

B  :  4-  <  KT  K)"1  KT(o-'  -  ar'Q)  (5  3) 

Unfortunately  the  true  stresses  involved  in  matrices  c'and  cr^  are  not  known.  A  satis¬ 
factory  approximation  to  B  can  be  obtained  as  follows: 
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In  other  words  the  coefficients  of  the  linear  stress  field  contained  in  B  are  obtained 
from  the  nominal  stresses  in  the  elements  in  the  region  of  point  0.  An  approximation  to  the 
true  stress  distribution  at  0  is  then  obtained  by  substituting  B  from  Equation  54  into 
Equation  45. 

crQ  =  a  ■+-  R  { <7  -  o-J )  (56) 


where  , 

R  =  -  -5-  5L  N 
0  B 

N  =  (KT  K)’1  KT 


(57) 

(58) 


EQUILATERAL  TRIANGULAR  GRID 


Figure  21.  Contiguous  Triangular  Elements 

Figure  21  shows  four  elements  formingpartof  an  equilateral  triangular  grid.  Corrections 
to  the  calculated  nominal  stresses  in  element  “0”  are  desired  based  on  the  assumption  that 
the  nominal  stresses  at  points  0,  1,  2  and  3  can  be  approximated  by  a  linear  stress  field. 
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Substituting  the  coordinates  of  the  four  points  shown  in  the  figure  into  the  Equation  for  K 
(Equation  50)  gives: 


K  = 


i  i  I  yr 


kt  k 


3 

0 

0 

-  3 
0 
0 

0 

0 

I  8 
0 

0 

0 


1 

1 

fL 

3 

1 

1  ° 

I  0 

1 

1 

1 

—  1 

1  -S3 

i 

1  0 

l 

4 

1 

i 

0 

r-v- 

[yr 

I 

i 

/3 

1  0 

1  0 

•V| 

3 


I  0  "|-2  /T 


-2/3  I 


I  0 

I  2/T 


0 

4 

0 

0 


0 

0 

4 

0 


0 

0 

0 

0 

0 

18 


(59) 


(60) 


The  matrix  K  K  is  not  only  nonsingular,  but  also  diagonal, 

T 

Substituting  K  and  K  K  into  Equation  58  gives: 


N 


0 

l 

6 

0 

1 

I 

0 

-1 

6 

0 

0 

0 

0 

0 

S* 

-1 

1 

1 

o 

S3 

, 

o 

-S3 

o 

12 

4 

1 

12 

4 

6 

1 

n 

-S3 

1 

1 

-1 

n 

-73 

n 

73 

S 

V 

12 

1 

1 

~4 

u 

12 

u 

6 

/3 

18 
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0 

1 

1 

i 

73 

18 

0 

0 

-VT 

9 

0 

0 

(61) 
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Equation  39  gives  the  nominal  uniform  stresses  in  the  triangular  element  designated  “0” 
in  Figure  21.  These  stresses  are  given  explicitly  by  Equations  26  for  a  triangle  located  as 
shown  in  Figure  17.  In  order  to  make  the  coordinate  axes  of  Figure  21  pass  through  the  centroid 
of  the  triangle,  take  =  -a/2,  y^  =  -  y3a/6.  Equations  26  then  become 


(62) 


Comparing  Equations  39  and  62  shows  that 


0  1 

73 

3 

0  1-73 

I 

1 

0  1 

273 

3 

! 

0  !  o 

V3  | 

-73  1 

2  1 

0 

~\  °  J 

Substituting  N  and  ov,  from  Equations  61  and  63  into  Equations  57  gives 

-D 


f 

t 

73 

1  1 

t 

73”  i 

1 

JL 

n 

6 

12 

1  2 

!  6 

”  12 

-,z  j 
-73"  1 

"  3 
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73 

i 

*  0 

_l_ 
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“  6 
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6 

6  ! 

3 

7L 

TIL 

-1 

iiTI 

-71 
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0 

0 

1 

24 

12 

24 

]  24 

12 

12 

(63) 


(64) 
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APPENDIX  HI 

EVALUATION  OF  THE  RECTANGULAR  BAR-PANEL  IDEALIZATION 


Figure  22  shows  the  coordinate  system  for  the  idealized  plate.  Figure  23  shows  the  notation 
for  displacements  and  strains  of  a  typical  set  of  bars  and  panels.  Denote  the  joint  whose 
coordinates  are  x.,  y^  as  the  ijth  joint.  In  the  following  definitions  the  expressions  “horizontal” 
and  “vertical,”  “left”  and  “right,”  “above”  and  “below,”  refer  to  directions  and  relative 
positions  of  bars,  panels  and  joints  appearing  in  Figure  23. 

Let 


u.,,  v., 
ij  ij 


=  components  of  the  displacement  of  the  ijth  joint. 


e 


=  the  tensile  strain  of  the  bar  joining  the  ijth  joint  and  the  shear  node  to  the  left 
of  this  joint. 


e 


/ 


=  the  tensile  strain  of  the  bar  joining  the  1-1,  jth  joint  and  the  shear  node  to  the 
right  of  this  joint. 


6y^  =  the  tensile  strain  of  the  bar  joining  the  ijth  joint  and  the  shear  node  below  this 

joint. 

./ 

€  y^  =  the  tensile  strain  of  the  bar  joining  the  i,  j-lst  joint  and  the  shear  node  above 
this  joint. 


=  the  shear  strain  of  the  panel  which  has  as  three  of  its  vertices  the  ijth  joint, 
the  i-1,  jth  joint  and  the  i,  j-lst  joint. 


a.,  b. 
i  J 


'll 


=  the  lengths  of  the  horizontal  and  vertical  edges  of  this  panel. 

=  the  horizontal  displacement  of  the  shear  node  to  the  left  of  the  ijth  joint. 


=  the  vertical  displacement  of  the  shear  node  below  the  ijth  joint. 
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DIFFERENCE  EQUATIONS 

Following  are  the  gpverning  equations  for  the  idealized  plate: 
STRAIN  -  DISPLACEMENT  EQUATIONS 
From  Figure  23, 


e‘ii 

=  ~  (u..  —  71..  ) 

0}  Ij  /|j 

(65) 

/ 

*x- 

XIJ 

■  -Vi,j> 

(66) 

€y.. 

~bj  (,'i  ”  tli  1 

(67) 

€iii 

-%%  -..1-1 » 

(68) 

Xxy.. 

7iJ 

) 

(69) 

STRAIN  RELATIONSHIPS 

E‘‘u 

=  c r - - - fb-  a 

xij  bj  +  b}+1  ^  J  yij 

+  bj+1 

(70) 

/ 

E£x 

_/  V  (  _ 

=  cr - lb. o' 

xi*H,j  bj  +bj+T  '  J  yi  j 

%.iJ 

(71) 

Ec 

yij 

—  v  /  — 

-  cr  -  -  (  o.  cr 

yij  ai  +ai+T  '  1  xij 

+  °j  +  l 

xi+l,j/ 

(72) 

/ 

-a  (a . c t 

yi,j+l  ai  +ai+1V  1  xij 

4-  a 

i+l 

a'  ) 

xi+l,  )' 

(73) 

Ey 

'xy., 

ij 

=  2(  14-  u)t 

xy.  • 

•  j 

(74) 
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JOINT  EQUILIBRIUM 


2~  ,0i  +°H-I  )t<Ty[  j+I 

1. 

?«•)  +  Vl"%  <?  — + 


2  (°i  +  0  m  ),<7y-j 


Figure  24.  Joint  Equilibrium 


From  Figure  24 


i+T,j  ij 


~ax.  =  0 


SHEAR  NODE  EQUILIBRIUM 


cr  -or  =  0 
yiJ+1  'j 


0 i f  r*y.  ... 
1 1,  J+l 


2(bj+bi+11,<,  ' 0— ►  +6  lHtr 


'J  * 


a .t  r 
i  xy . . 

'J 

Figure  25.  Shear  Node  Equilibrium 


Figure  25  shows  the  shear  node  to  the  left  of  the  ijth  joint. 


(75) 

(76) 


(77) 


67 


AFFDL-TR-68-150 


T  (°i  +ai+l)t<Tyfi 


b.tr 
j  *y: 


t 


u 


b.  iz 

J  xy. 


i+1  J 


2(ai  +ai+l)ta-y: 


Figuer  26.  Shear  Node  Equilibrium 


Figure  26  shows  the  shear  node  below  the  ijth  joint. 


— -{a .  +-  a 
2  1 


LINEAR  STRESS  DISTRIBUTION  BOUNDARY  CONDITIONS 


From  Equation  19,  neglecting  the  contribution  of  the  constant  stress  field: 


oj 


2b3x0+  6b4y 


°*x 


m  +!,  j 


2b3xm 


+  6b4yj 


y: 


10 


6tl  Xi  +2b2>0 


< 

y  i,  n  + 1 

=  6b.  xi  +2b2yn 

i 

Txyoj 

=  -2bzx0-2b3( 

T)C  y"i  +T,  j 

-  ~2bg  xm~  2b3  ^ 

v*> 

Txy. 

7I0 

1  ~  2  b  3  y0 

T*y. 

i,  n  + 1 

~  2Vn 

(78) 


(79) 


(80) 


(81) 
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SOLUTION,  LINEAR  STRESS  DISTRIBUTION 

The  following  displacements  and  stresses  satisfy  the  boundary  conditions  for  the  linear 
stress  distribution  and  the  difference  equations: 


Eu. . 
i  J 

=  -b  (3i/xz  4-  3y.z  )- 

i  i  'j 

Ev. . 
i  J 

=  b  (6x.  y.  )  +  b2(yj 

,  2 

Ei?..  =  -  b.  (3yx. 

i  i 

+  b3(xZ  - 

S  •  *>,<««,») 

v4 

2 


-b3(2i/x.  y.  -  i/b.  x.  )  -  b4(3vy*  +3x*  ~3i/b.y  .  ) 


2b3Xj  +6b4yj 
6b,  Xj  +■  2b2yj 


Txy..  =  -2b2(xi  -r)-2b3(yi-4-)  ) 


(82) 


(83) 


(84) 


The  displacements  and  stresses  given  by  equations  82  and  84  are  in  exact  agreement  with 
the  exact  displacements  and  stresses  given  by  equations  21  and  19,  when  the  results  derived 
from  the  structural  model  are  interpreted  as  follows: 


1.  Joint  displacements  are  equal  to  displacements  of  the  plate  at  corresponding  points. 

2.  Bar  normal  stresses  are  equal  to  plate  normal  stresses  at  points  in  the  plate  corresponding 
to  joints. 

3.  Panel  stresses  are  equal  to  plate  shear  stresses  at  points  in  the  plate  corresponding  to 
panel  centers. 

Shear  node  displacements  do  not  match  displacements  of  the  plate  at  corresponding  points. 
These  displacements  therefore  should  be  disregarded. 


69 


AFFDL-TR-68-150 


QUADRATIC  STRESS  DISTRIBUTION 


The  quadratic  stress  distribution 

crx  =  -  d,  (6x2  )  +-  d4(6xy  )  -  dg  (6x2  -12  y2) 

ay  =  -  d,  (6y2  -  I2x2)  +  d2(6xy)  -  ds  (6y2  ) 

Txy  =  d,  (12  xy  )- d2  (3x2)  -  d4(3y2)  +d5  (I2xy 

is  derived  from  the  stress  function 

<p  =  d)x4+d2x3y-3(d|+d5)x2y2  +  d4xy3+-d5y*  (8  6) 

The  stresses  given  by  Equations  85  satisfy  equilibrium  and  the  stress  function  given  by 
4 

Equation  86  satisfies  V  </>  =  0. 


The  plate  displacements  for  the  case 

d,  =  1  ,d2  =  d4  =  d5  =  0,  v  0  are 

E  u  =  -  2  x  3 
E  v  -  I  2  x2  y  -  2  y3 

The  displacements  for  the  case 

d2  =  I ,  d(  =d4  =  d5  =  0,  v  =0  are  v 

E u  s  “ y3 

Ev  =  3xy2  -  2x3  ' 


' 

► 

, 


BOUNDARY  CONDITIONS, 


d,  -- 

'»  d  2  =  d4  =  d5 

cr 

xoj 

=  -  6xZ 

0 

o’ 

*  .2,,  -6,J 

T 

Xyio 

1  12  1  «i  -  f-1  *0 

Txyoj 

■  ,2*„<  -  i> 

Similar  conditions  apply  to  the  other  boundaries. 


SOLUTION,  d,  =  1,  d2  =  d4  =  d5  =  0. 


(87) 


(88) 


(89) 
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The  following  displacements  and  stresses  satisfy  the  boundary  conditions  and  the  difference 
equations : 


Eu. . 

»J 

=  -2x? 

i 

—  a2x . 

E'u 

=  12  x  2 

y-  —  2  y  3 

J  J 

E^j 

=  -2x3 

i 

+  2  a  x  2 

i 

E^ 

2 

=  12  x 

i 

1  -  *>J 

-6bx2  +  3by2 
'  J 

% 

=  -  6xZ 

i 

°V.. 

7'J 

2 

=  1  2x. 

-  6yi 

T 
xy .. 

1 1 

=  12  (  x 

--b  ) 

2 

(90) 

(  91) 

(92) 

(93) 

(94) 


The  errors  obtained  by  comparing  these  solutions  with  the  exact  equations  are  as  follows: 


EAu, .  =  -  o2  x 

i  J 

EAv.. 

•J 


-  az  x .  ) 

! 


Acrx..  =  Aoy..  =ATxy.. 

ij  ij  Jij 


BOUNDARY  CONDITIONS,  dg  =  1,  d:  =  d4  =  dg  =  0 


=  0 


>oj 


=  6xi*o 

_ 

_  2 

T 

=  -  3x . 

xy. 

10 

r" 

■x  2 

=  —  3x . 

Xyoj 

\ 

(95) 

(96) 


(97) 


The  shear  stress  rxyj0  represents  the  average  shear  stress  at  the  shear  node  on  the 

lower  edge  of  the  structural  model  between  the  joint  at  x.  ..  n  and  the  joint  at  x.  .  The  shear 

2  l^-LjU  io 

stress  ^xy0j  equals  a  contribution  (-3xo  )  representing  the  constant  shear  stress  on  the  left 
hand  edge  of  the  model  and  a  contribution  (-a^)  necessary  to  provide  equilibrium  of  moments 
on  the  model.  Similar  conditions  apply  to  the  other  boundaries. 
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SOLUTION,  d_  =  1,  cL  =  d .  =  dc  =  0 
2  14  5 


The  following  displacements  and  stresses  satisfy  the  boundary  conditions  and  the  difference 
equations : 


Eu.  . 
i  J 

= 

E'u 

: 

,  2 

3  x .  y.  - 
'  J 

-  2x’  +-  b2 

X  , 

r 

P- 

UJ 

= 

3 

-yi 

ESU 

- 

-b  )-2x3 

i 

.2.  2 
+b  x. 

i 

% 

= 

0 

) 

% 

= 

6Vj 

1 

1 

► 

.. 

1 1 

-3  (xj  - 

□  f  a2 

J 

The  errors  are  as  follows: 

EAu  =  0  \ 

EAv  =  b2  x?  j 

i 

A°x  =  Acry  --  0  \ 


(98) 


(99) 


(100) 


(101) 


(102) 
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APPENDIX  IX 

STEADY  STATE  VIBRATION  OF  DAMPED  ELASTIC  STRUCTURES 


Equations  are  presented  for  the  analysis  of  steady  state  vibrations  of  damped  elastic 
structures  by  the  Matrix  Displacement,  Matrix  Force,  and  Unified  methods. 

GOVERNING  EQUATIONS 


Dynamic  Equilibrium: 

p_  Fr.  +  F„  +  P  F  +P„R+P 'id>  -  0  (103) 

E  E  D  D  mm  R  <p~ 

where  P£,  PD>  Pm>  P  R  and  P<fi  are  rectangular  matrices  of  components  in  the  joint 
degrees  of  freedom  of  unit  values  of  elastic  element  forces,  damping  element  forces,  inertia 
force  components,  reactions  and  external  loads,  and  Fg,  F^,  F  m>  R  and  <f>  are 
column  matrices  of  elastic  element  forces, damping  element  forces,  inertia  forces,  reactions 
and  external  loads.  The  later  group  of  matrices  are  time  dependent.  This  equation  can  also  be 
written  as  follows: 

PUEFE+PUDF0  +  PUmFm+PURR  +  V^  =  °  ('041 


PCEFE  +PCDFD+PCmFm+PCRR+V*>’  0  (I05J 

where  PyE’  ^UD'  ^  Urn*  ^UR  ant^  ^u4>  are  components  in  the  unconstrained  joint 
degrees  of  freedom  of  unit  values  of  the  elastic  element  forces,  the  damping  element  forces, 
the  inertia  forces  and  the  external  loads,  and  P^,  ^CD*  ^Cm’  ^CRanc*  ^c4>  are 
corresponding  matrices  for  the  constrained  joint  degrees  of  freedom.  The  unconstrained  joint 
degrees  of  freedom  are  selected  in  such  a  manner  that  reactions  have  no  components  in  these 
degrees  of  freedom. 

PUR  =  0  (106) 


FORCE-DEFORMATION  EQUATIONS 


ke, 


=  c  « , 


m 


mA 


m 


(107) 

(108) 
(  109) 
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where  k  ,  c  and  m  are  square  matrices  of  element  stiffnesses,  damping  constants  and 
masses  and  e  _,  e  „  and  A  _  are  column  matrices  of  elastic  element  deformations, 

hi  u  hi 

damping  element  deformations  and  joint  displacement  components,  respectively, 

EQUATIONS  OF  COMPATIBILITY 

The  following  equations  can  be  derived  from  the  principle  of  virtual  work: 


®E 

+  ®A 

=  -puea-4ar 

(  110) 

®D 

+  eB 

=  -PUD  A  -  Pcd  Ar 

(III) 

+ 

=  _PUm  A-pCm  Ar 

(112) 

where  A  and  A  R  are  column  matrices  of  components  of  joint  displacements  in  the  un¬ 
constrained  and  constrained  degrees  of  freedom  respectively  and  e  A  and  e  n  are  column 
matrices  of  unassembled  deformations  of  the  elastic  and  damping  elements  respectively.  The 
matrix  e  is  defined  subsequently. 

SINUSOIDAL  EXCITATION,  STEADY  STATE  RESPONSE 

Equations  104  to  112  inclusive  are  the  governing  equations  for  the  structure  subjected  to 
arbitrary  time  varying  external  loads,  support  displacements  and  unassembled  deformations. 
Assume  that  these  forcing  functions  are  sinusoidal  functions  of  time  of  a  single  frequency  u> 
but  varying  phase.  Any  such  function  Y  can  be  written  in  the  form 

Y  =  Yg  sin  ajt  +  Yc  cos  uit  (113) 

where  Yg  =  j  Y  |  cos  Q  Yc  =  j  Y  sin  6y  Y  j  is  the  amplitude  of  the  forcing 
function  and  6  is  the  phase  angle.  A  sufficiently  long  time  after  the  initiation  of  the  dynamic 
process  all  structural  responses  are  sinusoidal  functions  of  the  time  of  frequency  o>  .  Each 
response  therefore  can  be  written  in  the  form  of  equation  113.  In  the  following  the  quantities 
Y  and  Y  are  called  the  sine  and  cosine  components  of  Y.  Sine  and  cosine  components  are 
denoted  by  adding  subscripts  S  and  T  respectively  to  the  symbol  for  a  forcing  function  or  re¬ 
sponse.  For  example,  the  matrix  of  unconstrained  displacement  components  is  written 

A  =  AssinaU+Accosajt  (114) 


74 


AFFDL-TR- 68-150 


UNIFIED  MATRIX  METHOD 


Substitute  expressions  of  the  form  of  equation  114  for  F  F~,  F  and  d>  into 

ii  jj  m  * 

equation  104.  Equating  like  terms  yields 


where 


pu  Fs  +  PU^S  '  0 
pu  pc  +  pud>*c  =  0 


PU  =  .  PUE  PUm  PUd] 


FS  ={FES  Fm  S  FOS  } 


Fc  Tec  Fmc  Foc. 


and  the  symbols  £  j  and  |  |  denote  rectangular  and  column  matrices,  respectively. 


In  a  similar  manner  Equations  107,  108  and  109  yield 


Fo  =  C< 


where 


Fc  :Jes+  7)9. 


z 

■w  m 


=  {*ESA 


mS  #DSi 


{*  EC  ^ mC®  DC  } 


(120) 


Similarly  Equations  110,  111  and  112  can  be  written 

*S+  ®TS  s  “  PU  AS“  PC  ARS 

#C  +  #TC  =  “PU  AC-PCARC 


(121  ) 


where 


*TS  =  {  *AS  ®mS  * Bsl 
PC  =  [  PCE  PCm  PCd] 


#AC  *mC  eBCJ 


(122) 


(  123) 
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Let 


*■  . 


s 

'  f'  1 

c*  ■ 

0 

m 

z 

0 

CO 

c-' 

0 

CO 

:  i24) 


(125) 


Note  that  17*17  +  £  £  =  I  =  Z* V  *  0 

Solving  equation  118  for  eg  and  e^  and  substituting  the  results  into  Equation  121  yields 


7?*f54C‘fc+pJ  as+.ts+pJ  &rs=  0  ' 


s 

.* 


i  FS  +  ,j'FC  +  ,U  Ac+,rc  +fC  ARC  =  0 


(126) 


Equations  115  and  126  can  be  written  as  follows: 


’  1*  £*"  ■£  0  ' 

FS 

+ 

•is  +pc  ars 

-C*  1*  0  rTu 

Fc 

•tc  +pc  arc 

Py  0  0  0 

As 

pu <fc  0s 

1 - 

O 

TJ 

C 

O 

O 

1 _ 

1  > 

O 

*  - 

pu<£0c 

=  0 


(127) 


Equation  127  is  the  Unified  Matrix  equation  for  steady  state  damped  structural  vibration.  The 
equation  can  be  solved  for  the  sine  and  cosine  components  of  the  element  forces  and  dis¬ 
placements. 


ELASTIC  ANALOGY 


The  structure  can  be  replaced  by  an  analogous  structure  in  which  the  masses  are  replaced 
by  negative  springs.  Equations  119  show  that  the  unassembled  stiffness  matrix  for  these 
negative  springs  is  -  to  2m. 

Equation  120  shows  that  the  negative  spring  deformations  are  equal  to  the  joint  displace¬ 
ments.  The  negative  springs  therefore  must  be  considered  to  be  attached  to  the  frame  of 
reference.  Equations  122  show  that  ®mg  and  emQ  are  sine  and  cosine  components  of  un¬ 
assembled  deformations  of  the  negative  springs. 
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MATRIX  DISPLACEMENT  METHOD 


Solving  Equations  121  for  eg  and  e„  and  substituting  the  results  into  Equations  118 
yields 

F,  =  -  PS  A+  £  PI  Ac-  V'rs +£vc  - V  PJ  As+  £pcT  A 


FC  “  £PJAs  1)I>UAC~£,TS  V‘tc  £pc  As“,?PJ  A 

Substituting  Fg  and  Fc  from  Equations  128  into  Equations  115  yields 


RC 

RC 


S  — T  ~ 

AS 

= 

— *  — 

T  S 

AC 

fc 

where 


*s  5  Vs-  pu  ’Vs-%  £*tc  ~suc  As+trc  V 
%  Vc  -  pu  £*ts"pu  Vtc'  V  As-  *ucarc 


s  =  pu’?pu 


*pu£n 


u 


S  :  P 

uc  u 


*?pc  TucipuCpc 


(128) 


(129) 


(130) 


031) 


Equation  129  can  be  solved  for  the  sine  and  cosine  components  of  displacements  by  three 
methods.  Two  of  these  methods  involve  finding  S  “1  or  y  -1  or  solving  the  equations  by 
equivalent  procedures.  The  method  which  involves  S  ^  is  likely  to  be  inaccurate  in  the  case 
of  a  complicated  structure  having  many  natural  frequencies  because  S  is  poorly  conditioned 
when  the  impressed  frequency  is  close  to  a  natural  frequency.  The  matrix  T  is  likely  to  be 
poorly  conditioned  or  singular  if  some  of  the  structural  elements  have  zero  damping.  The 
third  and  best  procedure  involves  solving  the  complete  set  of  equations  by  elimination, 
employing  a  pivot  selection  technique.  Element  forces  can  be  calculated  from  Equations  128. 


Equations  128  and  129  are  the  Matrix  Displacement  equations  for  damped  structural  vibration. 
MATRIX  FORCE  METHOD 


The  restriction  that  reactions  have  no  components  in  the  unconstrained  degrees  of  freedom 
is  removed.  Equation  103  can  be  written 

pfrfr  +  Vs°  11321 
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where 

PRR  =  [PE  Pm  PD  Pr]  PR  ‘{FE  F»  PD  R  } 

The  element  forces  can  be  expressed  in  terms  of  the  redundants  and  the  external  loads  as 

Fr  :fXX+f^  C34) 

where  f  x  and  f  are  rectangular  matrices  of  element  forces  resulting  from  unit  values 
of  the  redundants  and  of  the  external  loads  respectively,  and  X  is  a  column  matrix  of  re¬ 
dundants.  Note  that  the  negative  springs  representing  the  masses  are  considered  part  of  the 
structure.  The  forces  in  these  negative  springs  can  be  selected  as  redundants.  The  matrices 
f  x  and  t<p  are  independent  of  time;  they  can  be  computed  from  P ^  and  P<£  through 
procedures  applicable  to  static  analysis. 

The  sine  and  cosine  components  of  F  are  obtained  by  substituting  expressions  of  the 
form  of  Equation  114  for  F  ,  X  and  <j>  into  Equation  134  to  yield 

frs  ='xV  V*s 
frc  •  fXXC  +  V^C 

The  compatibility  conditions  (Equations  126)  can  be  written 

*-  F  +  r-*  p  +pTA  +  * 

TR  rs  ^fr  rc  fr  “ucs  frs 

~  £  FR  FRS  ¥  ^FR  FRC  PFR  AUCC+  #  FRC 

where 

*-  r  *  “i  i.*  .  r  »•* 

17  FR  =  ^FR  ‘  ^ 

0  0 

^UCS  =  {^s  ars1  ^ucc  -  {  ac  arc1 

*  FR  S  {  eTS~  A  rS  }  ®FRC  =  {  *TC~ArC  } 

and  A  is  a  matrix  of  reaction  displacements. 


(137) 

(138) 
(  139) 
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Substituting  F  Rg  and  F  from  equations  135  into  equations  136  and  multiplying 
both  equations  through  by  f  T  gives 

X 


8xx  xs  +^xxxn+Sv 


where 


8XX  =  fX  %RfX 

\4>  z  f*  ^FR  V 

The  matrix  PFR  f  can  be  shown  to  be  null.  Therefore 


T  T 

R  fX}  AUCS  +  ,X  *FRS  =  0 

(140) 

frV  Aucc  +  ,x  *frcs0 

(141) 

I  T 

°XX  *  fX  *FR  fX 

(142) 

\<f>  :  fx  ^FR  f<£ 

(143) 

X 


^XX 

^XX 

1 

<s> 

X 

- 

r 

_  ^XX 

8XX^. 

1 

X 

0 

1 _ 

— 1 

0 

*< 

_ 1 

(144) 


where 


XS  =  8x^$  8X</>^c“  fX  eFRS 

Xc  =8X^S-Vc-fTX-FRC 


(14  5) 


Equations  144  can  be  solved  for  the  sine  and  cosine  components  of  the  redundants.  As  in 
the  Displacement  Method  the  best  approach  appears  to  be  to  solve  the  complete  set  of  equations 
by  elimination,  employing  a  pivot  selection  technique.  Element  forces  can  be  calculated 
from  Equation  135.  Deflections  can  be  obtained  by  calculating  the  deformations  of  the  negative 
springs  as  follows: 


AS  (jj2  PUm  m  mS  *Um  #mS 


kC  =  Je  PUmm  FmC  PUm*mC 


(146) 


Equations  135,  144,  and  146  are  the  Matrix  Force  equations  for  steady  state  damped 
structural  vibration. 
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